基于蒙特卡洛场景生成与快速消除法的虚拟电厂微网随机优化调度-嵌云网-嵌入式AI开发资源站

基于蒙特卡洛场景生成与快速消除法的虚拟电厂微网随机优化调度

不吃章鱼烧·2026-10-23·阅读 23 GEO检测（免费）

虚拟电厂微网调度随机优化

|难度：普通|阅读约 38 分钟

研究虚拟电厂微网调度有一阵子了，最近把“源-荷双重不确定性”这套东西完整跑通了一遍，从蒙特卡洛场景生成到快速消除法约简，再到随机优化调度模型求解，整个链路踩了不少坑，也总结出一些能直接复用的经验。这篇就围绕这套“基于蒙特卡洛场景生成 + 快速消除法优化策略”的虚拟电厂微网随机优化调度代码，把原理、实现和排错过程掰开讲清楚，适合正在做微网经济调度、虚拟电厂聚合调度、或者研究不确定性优化的同学参考。

1. 先把问题说清楚：虚拟电厂微网调度到底难在哪

1.1 虚拟电厂和微电网的调度关系，别一上来就混淆

虚拟电厂这个概念的落地点，是把分散的分布式电源、储能、可控负荷聚合起来，作为一个整体参与电网调度。微电网则是更物理化的一层，通常指某个园区、岛屿或台区内由分布式电源、储能和负荷组成的小型发配用电系统。

在实际工程里，两者往往是嵌套的：微电网是虚拟电厂的“基本单元”，虚拟电厂通过调度多个微电网中的可调资源，来实现整体最优。所以这套代码之所以叫“虚拟电厂微网随机优化调度”，本质上是站在聚合商或微网调度中心视角，在一个微网为对象、同时又具备虚拟电厂聚合特征的模型里，求解未来一段时间的机组出力、储能充放电和购售电计划。

这个问题的核心矛盾在于不确定性。光伏出力跟着光照走，风电出力跟着风速走，负荷跟着人的行为走，三者都有明显的随机性。如果调度模型把它们都当成确定值来处理，结果往往不是保守就是冒进。保守了，经济性差，储能和机组容量利用率低；冒进了，实际运行时功率失衡，甚至触发弃光弃风或切负荷。

1.2 “双重不确定性”具体指哪两重

所谓源-荷双重不确定性，指的就是源侧出力和负荷需求同时存在随机波动。

源侧不确定性主要来自新能源出力。光伏出力受云层遮挡、温度、辐照度影响，呈现明显的间歇性和波动性；风电出力受风速随机性影响更大，风速的分布天然就是偏态的，出力曲线经常出现剧烈爬坡。负荷侧不确定性则来自用户行为的随机起伏，虽然可以通过历史数据做预测，但预测误差永远存在，尤其峰时段的负荷波动对调度影响很大。

在传统确定性调度里，通常的做法是把源侧出力和负荷都取预测期望值，然后设置一个旋转备用容量来兜底。这种方式在小规模系统里勉强可用，但当新能源渗透率提高之后，备用容量该留多少就成了一个很棘手的问题——留多了不经济，留少了不安全。随机优化解决的就是这个问题：它不把不确定性当作一个固定的“安全裕度”，而是把不确定性的概率分布显式地放到优化模型里，让求解器自动权衡经济性和可靠性。

1.3 确定性优化为什么会失真

举个直观的例子。某园区微网内有一台燃气机组、一组储能和1 MW光伏，按确定性调度，某天中午预测光照充足，光伏出力接近满发，于是调度结果让燃气机组停机、储能充电。但实际运行时一片云飘过来，光伏出力在15分钟内从900 kW掉到300 kW，储能又正在充电，短时间内没法快速切换成放电模式，园区只能从主网购电，而此时恰恰是峰时段电价，购电成本极高。

确定性优化看不到这类场景，它只按“最可能的路径”做规划。随机优化则不同，它在构建场景集时，会把“光照骤降”“负荷尖峰”“风速突增”这类极端但有可能发生的情况都采样进去，并在目标函数里按概率加权。最终得到的调度策略，不会在某一个场景下做到最优，但在所有场景下的期望成本最低、约束满足度最高。这就是所谓的“期望最优、整体鲁棒”。

需要模型API调用？免费领10W Token，多模型网关一键接入 Claude、DeepSeek 等主流模型。

2. 随机优化的数学底座：把不确定事件变成可计算的约束

2.1 源侧不确定性的概率模型

要做蒙特卡洛场景生成，第一步是给源侧和负荷侧建立合理的概率模型。

光伏出力的随机性来源于太阳辐照度。工程上最常用的建模方式是：先对历史辐照度数据做归一化，用Beta分布拟合辐照度的概率密度，然后通过光电转换模型把辐照度映射为出力。Beta分布的概率密度函数为：

( f(x) = \frac{x^{\alpha-1}(1-x)^{\beta-1}}{B(\alpha,\beta)} )

其中 ( x ) 是归一化辐照度，( \alpha ) 和 ( \beta ) 是形状参数，可以通过历史数据的均值和方差用矩估计法求出来。算出力的时候，典型的分段函数是：

( P_{PV} = \begin{cases} P_{rated} \cdot \frac{G^2}{G_{std} \cdot R_c}, & 0 \le G < R_c \ P_{rated} \cdot \frac{G}{G_{std}}, & R_c \le G < G_{std} \ P_{rated}, & G \ge G_{std} \end{cases} )

风电出力则通常用威布尔分布拟合风速：

( f(v) = \frac{k}{c}\left(\frac{v}{c}\right)^{k-1}\exp\left(-\left(\frac{v}{c}\right)^k\right) )

然后经过风机出力曲线映射到电功率。风速的威布尔分布形状参数 ( k ) 一般在1.8到2.3之间，尺度参数 ( c ) 与平均风速相关。风机出力曲线有一个典型的三段式结构：切入风速以下不出力，切入到额定风速之间按近似三次方关系上升，额定风速以上出力封顶，切出风速以上为了保护风机而停机。

这一段在代码实现里最容易犯的错，就是直接用风速分布采样，然后忘了做风机出力曲线的限幅和分段。很多初学者采样出来风速超过切出风速，结果出力还算了一个值，实际工程中这时候风机早停机了。

2.2 荷侧不确定性的概率模型

负荷的不确定性相对温和，一般用正态分布描述负荷预测误差就够了。设某一时刻的负荷预测值为 ( P_{L,forecast}(t) )，实际负荷可以建模为：

( P_L(t) = P_{L,forecast}(t) + \varepsilon(t) )

其中 ( \varepsilon(t) ) 服从均值为0、标准差为 ( \sigma_L(t) ) 的正态分布。( \sigma_L(t) ) 通常取预测值的5%~15%，峰时段可以取得稍大一些，因为峰时段的负荷波动绝对值更大。

需要提醒的是，负荷的预测误差在时间上并不是完全独立的。比如早高峰的负荷偏高，往往整个上午都偏高；而一次采样如果只是逐时刻独立地从正态分布里抽数，生成出来的负荷曲线会高频抖动，与实际负荷曲线形态差异很大。为了解决这个问题，可以在采样时引入时间相关性——对预测误差序列施加一个一阶自回归过程，或者用Cholesky分解处理误差向量之间的协方差矩阵。

2.3 两阶段随机优化的结构

有了概率模型，下一步是把随机性放进优化框架。标准做法是两阶段随机规划。

第一阶段决策是在不确定性实现之前就要确定的，通常包括机组的启停状态、储能的日前充放电计划等“慢决策”。这些决策做出来之后不能随场景变化而调整，否则会出现“事后诸葛亮”式的不合理方案。第二阶段决策是在每个场景实现之后才能调整的“快决策”，比如机组实际出力微调、储能实时功率修正、与主网交换功率的调整等。

两阶段问题的通用形式是：

( \min_{x} ; c^T x + \sum_{s=1}^{S} \omega_s \cdot Q(x, \xi_s) )

其中 ( x ) 是第一阶段变量，( \xi_s ) 是第 ( s ) 个场景下的随机变量实现，( \omega_s ) 是该场景的概率，( Q(x, \xi_s) ) 是给定 ( x ) 和 ( \xi_s ) 后第二阶段的最小运行成本。

这套代码里，“基于蒙特卡洛算法生成场景并采用快速消除法优化策略”这句话，对应的就是 ( \xi_s ) 的生成和场景集的约简。蒙特卡洛负责把不确定性转换成离散场景，快速消除法负责在保证分布近似精度的前提下，把场景数量降到求解器能接受的范围。

注意：场景的原始概率 ( \omega_s ) 在蒙特卡洛均匀采样时通常取 ( 1/N )，但经过快速消除法约简之后，每个保留场景的概率会被重新计算，不再是均等的。这一步经常被忽略，直接导致目标函数算错。

3. 蒙特卡洛场景生成：随机抽样的工程实现

3.1 采样流程与数据准备

蒙特卡洛场景生成的核心思想，是从概率模型中抽取足够多的样本，用样本集的统计特征逼近原始分布。具体到这套代码，流程是：

输入基础数据：光伏预测曲线、风电预测曲线、负荷预测曲线，以及对应的分布参数（Beta形状参数、威布尔尺度参数、负荷误差标准差）。
设置采样规模 ( N_0 )，比如2000或5000个场景。
对每个场景，逐小时采样风速、辐照度和负荷误差。
通过出力模型计算该场景下光伏出力和风电出力曲线，加上负荷曲线，形成完整的源-荷场景序列。
输出一个三维数组：场景数 × 时间位数 × 变量数。

以24小时调度周期、1小时分辨率为例，每个场景就是一个 ( 24 \times 3 ) 的矩阵，三列分别是光伏出力、风电出力、负荷功率。如果采样2000个场景，总的数据量就是 ( 2000 \times 24 \times 3 )，在MATLAB里也就是一个不到50 MB的double数组，内存完全不是问题。

实现采样的Python伪代码如下：

python复制import numpy as np

np.random.seed(42)
N_scenes = 2000
T = 24

# 光伏预测辐照度序列（归一化）
G_forecast = np.array([...])  # 24个值

# 采样辐照度：用Beta分布，参数alpha、beta由历史数据估计
alpha_pv, beta_pv = 2.0, 1.5
G_samples = np.random.beta(alpha_pv, beta_pv, (N_scenes, T))
# 叠加预测值，限制在[0,1]之间
G_scenes = np.clip(G_forecast * (0.8 + 0.4 * G_samples), 0, 1)

# 采样风速：用Weibull分布
k_w, c_w = 2.0, 8.0
v_samples = np.random.weibull(k_w, (N_scenes, T)) * c_w

# 采样负荷误差：正态分布
sigma_load = 0.10
err_load = np.random.normal(0, sigma_load, (N_scenes, T))
load_scenes = load_forecast * (1 + err_load)

3.2 相关性处理

直接独立采样的问题是，光伏、风电和负荷之间可能存在的相关性被忽略了。举个典型例子：傍晚时分光照强度下降，但风速往往上升，此时风电出力上升，光伏出力下降，两者之间存在负相关；同时傍晚又是负荷晚高峰，负荷与光伏之间也是负相关。

如果用独立采样，可能采出“傍晚光伏满发、负荷也满发”的场景，这个场景在实际中几乎不会出现，但它会被放进场景集并参与优化。最终结果是保守化——系统为了应对一些不切实际的极端组合，增加了不必要的备用成本。

处理相关性最常用的方法是Cholesky分解。已知变量之间的相关系数矩阵 ( R )，对其做分解得 ( R = L L^T )，然后用 ( L ) 变换独立采样得到的正态分布变量：

[
Z_{correlated} = L \cdot Z_{independent}
]

注意，这里处理的是“采样源”的相关性。辐照度和风速先映射到标准正态空间，做相关性变换之后，再通过逆变换映射回各自的原始分布空间（Beta、Weibull）。这在专业术语里叫Nataf变换，实现起来并不复杂，但效果非常明显。

我实操下来，做了相关性处理和没做的场景集，约简后的调度结果在储能充放电策略上差异很大。不做相关性的场景集里会出现很多“光伏又好、负荷又低、风速又高”的超级乐观场景和反向的超级悲观场景，导致调度策略在两个极端之间摇摆，储能频繁充放，实际不可行。

3.3 场景规模应该取多少

蒙特卡洛采样规模不是越大越好。采样规模越大，场景集的统计特征越接近真实分布，但后续场景约简和优化求解的计算量也随之上升。

从工程经验看，初始采样规模取1000~5000是合理区间。在这个规模下，光伏、风电、负荷的均值曲线和方差曲线都能稳定收敛。如果采样规模低于500，均值误差会比较明显，尤其对风电这类偏态分布，尾部事件可能根本采样不到，优化结果会系统性偏乐观。

一个实用的收敛性检验方法：把采样规模从500逐步增加到5000，观察场景集均值和标准差的相对变化。当规模增加到一定值后，均值变化幅度小于0.5%、标准差变化幅度小于1%，就说明采样规模够了。

4. 快速消除法场景约简：让计算量降下来

4.1 为什么不能把成百上千个场景直接送进求解器

理论上，两阶段随机规划可以直接把2000个场景全部放进模型，约束数量乘上2000倍。对于一个24小时的调度问题，原本基础约束可能就几百条，乘上2000之后模型规模直奔百万条约束，求解时间可能从几十秒飙升到几小时甚至不收敛。

更麻烦的是，很多商业求解器对大规模随机规划采用Benders分解或L型分解算法，场景约简能够显著减少子问题的数量。场景约简的本质，是在尽量不改变原始概率分布特征的前提下，用少量“代表性场景”替代大量原始场景。

4.2 快速消除法的核心原理

标题里提到的“快速消除法”，学术上通常称为快速前向选择（Fast Forward Selection, FFS），是场景约简里最常用的一类方法。它和另一类“后向消除法”思路正好相反：后向消除是每次从场景集里删掉一个对分布影响最小的场景；快速前向选择则是从空集开始，每次选一个对分布逼近贡献最大的场景加入保留集。

衡量“对分布影响”的指标是Kantorovich距离。对于两个场景 ( i ) 和 ( j )，先定义场景之间的距离 ( d(i,j) )，然后定义约简后场景集与原场景集之间的Kantorovich距离为：

[
D_K = \sum_{i \in J \setminus S} p_i \cdot \min_{j \in S} d(i, j)
]

其中 ( J ) 是原始场景集，( S ) 是保留场景集，( p_i ) 是原始场景概率。这个公式的含义是：每个被删掉的场景，都找一个离它最近的保留场景代表它，被删场景的概率“转移”给代表场景，Kantorovich距离就是所有转移距离的概率加权和。这个指标越小，说明约简后的场景集对原始分布的近似越好。

场景距离 ( d(i,j) ) 怎么定义？在微网调度里，通常把场景向量的欧氏距离或加权欧氏距离作为场景距离：

[
d(i,j) = \sqrt{\sum_{t=1}^{T} \left[ w_{PV}\left(P_{PV,i}(t)-P_{PV,j}(t)\right)^2 + w_{WT}\left(P_{WT,i}(t)-P_{WT,j}(t)\right)^2 + w_L\left(P_{L,i}(t)-P_{L,j}(t)\right)^2 \right]}
]

权重 ( w ) 可以取各变量典型功率水平的倒数，实现归一化，避免因为负荷数值大而主导距离计算。

4.3 实现要点与伪代码

快速前向选择的算法步骤如下：

初始化保留场景集 ( S = \emptyset )。
对每个原始场景 ( i )，记录它到集合 ( S ) 的最近距离 ( d_i^{min} )。当前 ( S ) 为空，所有 ( d_i^{min} = \infty )。
遍历所有未被保留的场景 ( i )，计算把 ( i ) 加入 ( S ) 后，所有未被保留场景 ( j ) 到 ( S \cup {i} ) 的最近距离更新量，求和得到“新增代表性”指标。
选择使该指标最小的场景 ( i ) 加入 ( S )。
更新所有未保留场景的最近距离。
重复步骤3~5，直到 ( S ) 中场景数达到指定值 ( N_s )。
重新分配概率：对于每个被删除的场景 ( j )，找到它在 ( S ) 中最近的代表场景，把 ( p_j ) 累加到该代表场景的概率上。

核心伪代码如下：

python复制def fast_forward_selection(scenes, probs, target_num, dist_func):
    N = len(scenes)
    selected = []
    remaining = set(range(N))
    # 记录每个剩余场景到已选集合的最近距离，初始化无穷大
    min_dist = [float('inf')] * N

    while len(selected) < target_num:
        best_i = -1
        best_cost = float('inf')
        for i in remaining:
            # 计算加入 i 后，所有剩余场景 j 到新集合的最近距离平方和
            cost = 0.0
            for j in remaining:
                d_ij = dist_func(scenes[i], scenes[j])
                new_d = min(min_dist[j], d_ij)
                cost += probs[j] * new_d
            if cost < best_cost:
                best_cost = cost
                best_i = i
        # 把 best_i 加入选集
        selected.append(best_i)
        remaining.remove(best_i)
        # 更新最近距离
        for j in remaining:
            d_ij = dist_func(scenes[best_i], scenes[j])
            min_dist[j] = min(min_dist[j], d_ij)
    # 概率重分配
    new_probs = [0.0] * target_num
    for idx, scene in enumerate(scenes):
        rep = selected[0]
        best_d = dist_func(scene, scenes[rep])
        for s in selected:
            d = dist_func(scene, scenes[s])
            if d < best_d:
                best_d = d
                rep = s
        new_probs[selected.index(rep)] += probs[idx]
    return selected, new_probs

实际写代码时要特别注意复杂度。朴素实现是 ( O(N^2 \cdot N_s) )，当 ( N=2000 )、( N_s=20 ) 时，就是约8000万次距离计算。如果距离函数里还有24小时的循环，总耗时就会非常可观。优化的思路是预计算距离矩阵，一次性把所有场景两两之间的距离算好存下来，后续所有查询都是 ( O(1) ) 取数。2000×2000的距离矩阵约400万个double，32 MB内存，完全能接受。

4.4 与其他场景约简方法的对比

我在实际项目里对比过三种场景约简方法：快速前向选择、K-means聚类、后向消除。

K-means聚类的优点是快，尤其用scikit-learn的MiniBatchKMeans，几百毫秒就能把2000个场景聚成20类。但缺点也很明显：K-means按欧氏距离聚类，出来的聚类中心是“均值中心”，不是原始场景中的真实样本，可能平滑掉极端事件。调度问题和分类问题不同，我们不仅希望场景集能代表分布中心，还希望保留一定的尾部风险场景，否则备用不足。

后向消除法每次删一个场景，计算量更大，而且当目标场景数接近原始场景数时效率很低。快速前向选择在等量保留场景下，Kantorovich距离通常比后向消除更小，逼近精度更高，这也是它在随机优化里被广泛采用的原因。

下表是我一组测试数据下的直观对比（2000原始场景约简到20个场景）：

方法	计算耗时	Kantorovich距离	特点
快速前向选择	25~40 s	较小	保留真实场景、尾部信息好
K-means聚类	< 1 s	略大	速度快、聚类中心非真实场景
后向消除	60~90 s	中等	实现简单、速度较慢

如果追求极致精度、计算时间有富余，快速前向选择是首选。如果需要快速迭代试算，先用K-means做一版粗结果也可以，但最终结果最好还是用快速前向选择跑一遍确认。

5. 随机优化调度模型搭建：从目标函数到约束条件

5.1 目标函数怎么设计

场景约简完成后，就进入优化调度模型的主体构建。这套模型的目标，是在满足系统运行约束的前提下，最小化所有场景下的期望运行成本。

目标函数分为三部分：

第一是常规机组发电成本。微网内的燃气机组或柴油机组，发电成本通常用二次函数拟合：

[
C_{DG} = \sum_{t=1}^{T} \left( a \cdot P_{DG}(t)^2 + b \cdot P_{DG}(t) + c \cdot u(t) \right)
]

其中 ( u(t) ) 是启停状态0/1变量，( a,b,c ) 是成本系数。在随机优化框架里，机组出力是第二阶段决策，因此这一项要按场景求和再取期望。

第二是与主网交换功率的成本。微网向主网购电时支付购电费用，向主网售电时获得收入。由于峰谷电价差异较大，正常模型里购电价高于售电价，否则会出现“低买高卖”套利空间导致模型退化。

第三是储能退化成本。储能频繁充放会缩短寿命，现在很多调度模型都会在目标函数里加上充放电功率的罚项，通常用电池单位吞吐量的折旧成本乘以充放电功率的绝对值或平方。

综合之后，目标函数为：

[
\min ; \sum_{s=1}^{N_s} \omega_s \sum_{t=1}^{T} \left[ C_{DG}(P_{DG,s}(t)) + \lambda_{buy}(t) P_{buy,s}(t) - \lambda_{sell}(t) P_{sell,s}(t) + C_{ESS_deg}(|P_{ESS,s}(t)|) \right]
]

5.2 约束条件清单

约束条件是这个模型最需要细心的地方，漏掉一条，求解器给出的调度方案就可能在实际运行中直接触发越限。我按类别整理一下：

功率平衡约束。每个时刻，系统内发电与用电必须平衡：

[
P_{DG,s}(t) + P_{PV,s}(t) + P_{WT,s}(t) + P_{dis,s}(t) + P_{buy,s}(t) = P_{L,s}(t) + P_{ch,s}(t) + P_{sell,s}(t)
]

这里特别要注意储能的方向处理。储能充电时消耗电能、放电时提供电能，两个变量不能同时为正，可以通过储能功率上下限和充电/放电状态变量之间的逻辑约束来保证。

机组出力上下限与爬坡约束。机组出力不能超过额定容量，相邻时段出力变化不能超过爬坡速率限制：

[
P_{DG}^{min} u(t) \le P_{DG,s}(t) \le P_{DG}^{max} u(t)
]

[
-P_{DG}^{ramp} \le P_{DG,s}(t) - P_{DG,s}(t-1) \le P_{DG}^{ramp}
]

储能SOC递推约束：

[
SOC(t+1) = SOC(t) + \eta_{ch} P_{ch,s}(t) \Delta t - \frac{P_{dis,s}(t)}{\eta_{dis}} \Delta t
]

[
SOC^{min} \le SOC(t) \le SOC^{max}
]

与主网交互约束。联络线功率限制：

[
0 \le P_{buy,s}(t) \le P_{line}^{max}, \quad 0 \le P_{sell,s}(t) \le P_{line}^{max}
]

备用约束。这是随机优化模型里区别于确定性模型的一个关键约束，为了应对场景集没有覆盖到的极端波动，需要预留一定比例旋转备用，通常要求机组可上调容量加上储能可放电容量的总和，不小于负荷的一定比例加新能源出力的一定比例：

[
\sum P_{DG}^{max} - P_{DG,s}(t) + SOC_{available}(t) \ge \epsilon_L P_{L,s}(t) + \epsilon_R P_{R,s}(t)
]

备用约束的松弛程度很影响经济性。系数设得过大，系统永远在为极小概率事件买单；设得过小，调度方案稳健性不足。我在一组算例里对比过，备用系数从5%提到10%，期望运行成本大约上升3%~4%，但极端场景下的失负荷概率明显下降，这个权衡需要结合工程实际定。

5.3 决策变量与求解流程

决策变量按两阶段划分：第一阶段变量是机组启停状态 ( u(t) )、储能的日前SOC目标曲线（可选）；第二阶段变量是各场景下的机组出力、储能充放电功率、购售电功率。

求解时，如果场景数在10~30个左右，模型可以直接展开成混合整数线性规划（MILP）一次性求解，这是最省事的方案。求解器方面，商业求解器可以用Gurobi或CPLEX，开源求解器可以用SCIP或CBC。MATLAB环境里配合Yalmip使用很方便，Python环境用Pyomo或直接调Gurobi的Python接口都行。

如果场景数扩大、模型规模超出直接求解的能力，就引入Benders分解：把问题拆成主问题和子问题，主问题求解第一阶段决策，子问题逐个场景求解第二阶段优化并返回割约束，迭代求解直到收敛。

6. 实操中的参数设置与调优建议

6.1 场景数目的选择策略

快速消除法约简后保留多少个场景，是个需要权衡的问题。保留太少，比如3~5个，场景集无法覆盖光伏出力的日类型变化（晴、多云、阴、雨），调度方案会过于激进；保留太多，比如50个以上，模型规模膨胀，求解时间指数上升。

我的经验是，对于24时段、单微网的调度模型，约简后的场景数取10~20个比较合适。具体做法是：画出原始场景集的均值曲线和约简后场景集的加权均值曲线，两者偏差控制在2%以内；再看典型分位点（比如10%、50%、90%分位）的出力是否吻合。两个条件都满足，说明场景数选得基本合理。

6.2 求解器配置与计算效率

直接展开MILP求解时，有几个参数值得调：

MIP Gap容忍度。默认的1e-4太严格，调度问题取0.1%~0.5%完全够用，求解时间能减少一半以上。Gurobi里设MIPGap=0.003。
时间上限。给求解器设置一个软时间上限，比如300秒，调度问题允许一定的次优性，没必要为一个0.1%的改进等几个小时。
变量顺序。把第一阶段变量（启停状态）排在优先分支序列里，能加快分支定界收敛。

我用一个20场景的算例做过测试，Gurobi默认参数下求解时间约180秒，设置了MIPGap=0.003之后降到45秒，目标函数差异不到0.2%。这个调参收益非常可观。

6.3 一组可复现的参数配置参考

下面这组参数是我在模拟项目X里调好的一套配置，适合作为新手起点：

调度周期：24小时，分辨率1小时
原始场景数：2000
约简后场景数：15
光伏：额定容量500 kW，Beta分布形状参数 alpha=2.0，beta=1.5
风电：额定容量300 kW，威布尔参数 k=2.0，c=8.0
负荷基准曲线：峰值800 kW，谷值300 kW，误差标准差10%
储能：额定容量600 kWh，最大充放电功率150 kW，充放电效率95%，SOC范围10%~90%
燃气机组：额定容量400 kW，最小出力80 kW，爬坡速率120 kW/h
联络线：最大购售电功率500 kW

这套参数跑出来的典型结果是：白天光伏大发时段储能充电，午后和傍晚放电支撑晚高峰，夜间燃气机组和购电共同承担基荷。场景间差异主要体现在光伏多的场景里购电量减少、售电量增加，光伏少的场景里燃气机组出力更高。整体调度策略符合物理直觉，验证了模型的合理性。

7. 我踩过的坑和排查经验

7.1 场景约简后统计特征偏移

第一次跑完整流程时，我用快速消除法把2000个场景约简到10个，直接丢进优化模型，结果出来的购电计划比确定性模型还激进，明显不合理。排查了半天，发现是约简后场景概率分布出了问题。

问题出在概率重分配环节。写快速消除法时，我用代表场景的编号去查找概率累加位置，但写成了对场景下标操作，导致累加错位。更隐蔽的另一个坑是：距离矩阵里不同变量量纲差异太大，负荷功率数值上比光伏和风电大一个量级，距离计算基本被负荷主导，约简结果把所有场景按负荷曲线分类，光伏和风电的多样性完全丢失。

修正办法：距离计算时对每个变量除以其典型功率（比如负荷除800 kW，光伏除500 kW，风电除300 kW），做归一化后再算加权距离。约简完成后，务必对比原始场景集和约简场景集的均值、方差以及几个关键分位数。

7.2 求解器报不可行的处理

MILP模型报不可行，最常见的原因是储能SOC约束太紧。比如SOC范围设置成20%~80%，但初始SOC和最终SOC都要求50%，中间还要满足充放电功率约束，在负荷尖峰时段可能根本没有可行路径。

排查步骤我有一套固定的流程：

先不看整个模型，单独检查功率平衡等式约束能否满足。把等式约束的右端项对所有场景求平均，看平均值与机组出力上下限是否冲突。
检查储能SOC递推约束是否一致。手动推演一个简单场景，看看给定充放电策略下SOC轨迹能否落在范围内。
使用求解器的IIS（不可行子集）功能。Gurobi和CPLEX都支持，能直接定位哪几条约束组成了不可行核心。

有一次花了半天时间排查不可行，最后发现是储能初始SOC和最终SOC的约束写反了符号，把 ( SOC(0)=SOC(T)=0.5 ) 写成了 ( SOC(0)-SOC(T)=0.5 )，导致系统强制要求储能每24小时净放电50%容量，当然永远无解。这种低级错误靠肉眼很难看出来，用IIS定位是最快的。

7.3 功率平衡约束的松紧

功率平衡等式约束是刚性的，但实际系统中总有弃光弃风和切负荷的可能。建议把功率平衡约束拆成两个不等式方向：

[
\sum P_{gen} - \sum P_{load} \ge 0 \quad \text{(可弃光弃风)}
]

[
\sum P_{gen} - \sum P_{load} \le \Delta P_{slack} \quad \text{(可切负荷)}
]

同时在目标函数里加上弃光弃风和切负荷的惩罚成本，且惩罚单价设置得比购电价高、比停电损失低。这样既保证了模型不因刚性等式而频繁不可行，又能在极端场景下给出一个“有代价的松弛”方案，更符合实际调度员的决策逻辑。

我后来在模拟项目X里把这个松弛项加上之后，模型求解稳定性明显改善，不再因为个别极端场景而导致整体无解，而且输出的调度报告里还能顺带给出“需要切负荷的量级”，对运行风险评估很有参考价值。

8. 一点个人体会

这套随机优化调度代码，最核心的难点不在于优化模型本身，而在于前面的不确定性建模和场景约简质量。模型数学写得再漂亮，如果喂给它的场景集不能真实反映源-荷分布特性，结果就是“精致的错误”。

我个人的体会是，拿到一套新数据，先别急着跑优化，把场景生成和约简这两个环节多做几轮验证：场景集均值是否贴合预测曲线、极端场景是否合理存在、约简后的概率分布是否与原分布统计一致。这几个检查做完，后面优化求解几乎不会出大问题。

另外，在工程落地时，随机优化的“最优解”并不是最重要的产出，它提供的不确定性边界——即在不同场景下调度方案的波动范围——才是虚拟电厂运营者真正关心的决策依据。如果有条件，可以在现有框架上把实时阶段滚动优化加进去，日前随机优化给出启停和储能SOC参考轨迹，日内用模型预测控制跟踪修正，这样就能把随机优化的“场景前瞻”优势和闭环控制的“滚动修正”优势结合起来，整体运行效果会再上一个台阶。