1. 微网容量配置到底在解决什么问题
做微电网规划的人应该都清楚,容量配置这件事看着简单,真正落地的时候到处都是坑。风光机组装多少、储能配多大、柴油机和电网交互容量取什么值,这些决策直接决定项目全生命周期的经济性。但最麻烦的是,微网里的负荷和新能源出力都是不确定的,你按确定性优化算出来的配置方案,碰上极端天气或者负荷尖峰,可能直接扛不住。
两阶段鲁棒优化算法就是用来处理这种“先决策、后运行”问题的经典框架。我这里说的“两阶段”不是指两个求解步骤,而是指决策变量的两层结构:第一阶段做容量规划,第二阶段在给定容量下做最坏场景的运行调度。两阶段鲁棒优化的核心思想是:我不管不确定参数具体取什么值,我只保证在不确定集合内所有可能场景下,系统都能安全运行,并且运行成本不超过某个上限。这个思路在微网电源容量配置里特别合适,因为规划设计阶段你不可能预知未来每一刻的风速和光照,但你可以圈定一个合理的波动范围。
这篇代码做的事情,就是用MATLAB把整个建模、求解、迭代过程串起来,输入微网的负荷曲线、风光出力历史数据、设备参数和不确定集合设置,输出最优的电源配置容量和对应的运行成本。接下来我按自己的理解,把整个代码的设计思路、关键实现细节和调试过程中踩过的坑都拆开讲一遍。
2. 数学模型:容量配置问题的完整表述
2.1 目标函数与决策变量拆分
两阶段鲁棒优化在微网容量配置里的标准数学形式,可以写成这样一个紧凑的问题:
min C_inv(x) + max_{u∈U} min_{y∈F(x,u)} C_ope(y)这个式子看着简短,实际上分了三层结构:
- 外层min是第一阶段决策,变量x代表各类电源的安装容量,比如光伏的kW数、风机的台数、储能系统的kWh容量和kW功率、柴油机台数。C_inv(x)是这些设备的投资成本,按等年值系数折算到每一年的费用,这样才能和年运行成本放在同一个量纲下比较。
- 中间的max是鲁棒优化的核心,它在不确定集合U里寻找使运行成本最大的那个“最坏场景”,这个场景通常对应负荷高涨同时新能源出力低谷的情况。
- 内层min是第二阶段运行优化,变量y是每个时段的机组出力、储能充放电功率、弃风弃光量、切负荷量等,约束包含功率平衡、机组出力上下限、储能SOC递推等式、线路传输容量等。
把投资和运行分开处理,是因为两类决策的时间尺度完全不同。设备容量一旦定下来,二十年都不大变,但运行调度是每个小时甚至每分钟都要做的决策。你不可能用一个单层优化把所有时间尺度都揉在一起,那样问题规模和求解难度都会失控。
2.2 不确定集合的构建方式
鲁棒优化的结果是否合理,很大程度上取决于你怎么样描述不确定性。常用的集合是盒式不确定集合,也就是每个不确定参数独立地在区间内波动:
U = { u: u_min ≤ u ≤ u_max }具体到微网场景,光伏出力和负荷的不确定参数可以写成预测值加减波动偏差的形式。比如光照预测为100kW,允许偏差±20%,那么不确定集合就是[80, 120]kW。盒式集合的好处是建模简单、求解方便,但它的问题在于过于保守——它允许所有参数同时取到最坏值,而现实中这种情况几乎不会发生。
所以代码里通常会引入预算不确定集合,也就是限制所有不确定参数偏离预测值的总量不能超过一个预算值Γ。这个Γ的含义可以理解为“最坏情况出现的程度”,Γ越大结果越保守,Γ等于0时退化为确定性优化。理论上有名的Bertsimas和Sim方法就是干这个事的,它的好处是给了决策者一个调节保守度的旋钮,你可以根据项目风险偏好来选。
在MATLAB代码实现里,不确定集合的构建方式决定了后续生成极端场景的具体形式。我建议写成独立函数,输入是预测值序列和偏差比例,输出是一个描述盒式或预算集合的结构体。这样后续换保守度只需要改一个参数。
3. 两阶段鲁棒求解框架:列约束生成的核心逻辑
3.1 为什么要用CCG而不是Benders
两阶段鲁棒优化问题的求解方法,主流就两个:Benders分解和列约束生成。早期文献里Benders用得比较多,但它处理第二阶段是线性规划问题时需要引入对偶变量,然后在大M法处理双线性项时很容易出现数值问题。C&CG(Column-and-Constraint Generation)的好处是直接把第二阶段的场景作为新变量加入主问题,迭代过程中主问题规模会变大,但收敛速度明显更快,而且在处理混合整数第二阶段问题时不需要对偶转换。
这套代码选用的就是C&CG框架。整个流程可以分解为以下几个步骤:
- 初始化:设下界LB=-inf,上界UB=+inf,迭代次数k=1,初始最坏场景可以直接用预测值场景。
- 求解主问题:主问题的形式是在已知一组离散场景下做联合优化,目标函数是投资成本加上这组场景下运行成本的最大值,约束包含容量约束和每个场景对应的运行约束。主问题的最优值更新下界LB。
- 求解子问题:固定第一阶段变量x的取值,在不确定集合中寻找使运行成本最大的场景u,同时计算对应的运行成本。子问题的最优值加上投资成本更新上界UB。
- 收敛判断:如果(UB-LB)/UB小于设定间隙(比如1%),停止迭代;否则把子问题找到的最坏场景作为新一列加入主问题的场景集合,转回步骤2。
这个迭代过程里,最关键的环节是第3步怎么高效求解那个max-min双层问题。
3.2 子问题的对偶变换与KKT条件处理
子问题需要对每一个时段t和每一个不确定参数分别求解最坏场景,但直接写max-min是没法用现成求解器解的。标准做法是取内层min的对偶,把问题转化为单层max问题。
拿最简单的情况举例,假设第二阶段运行优化是线性规划,目标是最小化运行成本cTy,约束是Ay ≤ b + Bu,那么取对偶后可以写成:
max (b + Bu)^T λ s.t. A^T λ ≤ c λ ≥ 0λ是对偶变量。这时候目标函数里出现了u^T B^T λ,这是双线性项。好消息是不确定参数u是连续变量且约束通常是盒式的,所以这个双线性项和λ的乘积在u取区间端点时可以取得极值。具体来说,如果B^T λ的系数为正,u就取区间上界;如果系数为负,u就取区间下界。这就把max问题转成了一个混合整数线性规划,要么用大M法线性化,要么直接按系数符号判断,在代码里用逻辑索引赋值就行。
但要注意,第二阶段运行为线性规划这个前提很重要。如果你的模型包含储能SOC的0/1状态变量或者机组的启停变量,第二阶段就变成混合整数规划,对偶这条路就走不通了。这时候主流方案有两个:一个是把第二阶段建模成线性规划,储能和机组都用连续变量近似;另一个是改写成博弈论形式的KKT条件,把内层问题的最优性条件写成互补松弛约束,再引入大M法线性化。后者的处理复杂度会明显上升,但能保留整数变量的建模精度。
我写的代码里第二阶段用连续变量建模,储能SOC和充放电功率解耦成两个变量,通过SOC递推方程关联。这样保证子问题的线性规划性质,CCG迭代可以稳定收敛,单次子问题求解速度也快。
4. MATLAB代码实现与关键环节解析
4.1 代码的整体结构
这套代码按功能模块拆分为6个部分,运行时按顺序调用:
| 模块 | 文件/函数名 | 核心职责 |
|---|---|---|
| 数据输入 | load_case_data.m | 读取负荷、光伏、风机时序数据,设备参数表 |
| 主问题建模 | build_master_problem.m | 构建投资决策变量和场景相关的运行变量,返回YALMIP优化对象 |
| 子问题建模 | build_subproblem.m | 构建固定容量下的最坏场景寻优模型 |
| 不确定集合 | build_uncertainty_set.m | 生成盒式/预算不确定集合,提供场景生成接口 |
| 求解主程序 | main_ccg_solve.m | 执行C&CG迭代,输出配置方案和迭代收敛曲线 |
| 后处理 | post_process.m | 绘制容量配置图、运行曲线、迭代间隙变化 |
在主循环里,主问题用YALMIP建模,求解器调cplex或者gurobi都可以。子问题稍微麻烦点,因为它要显式操作对偶变量,YALMIP的dual函数可以自动提取对偶信息,但如果你想要更透明的控制,可以直接用YALMIP的灵巧建模写出max形式然后调用求解器。
4.2 主问题的建模细节
主问题的目标函数是投资成本,约束是容量上下限和一组已知场景下的运行可行性。每一次CCG迭代后,主问题里会多一个场景块,也就是多一组运行变量和运行约束。这里有一个容易踩坑的地方:运行变量必须显式区分场景索引,不能在不同场景之间共享第二阶段的决策变量,否则会破坏“最坏场景下可运行”的逻辑。
举个例子,储能SOC变量要写成soc(i, t, k),i是节点编号,t是时段,k是场景编号。不同k之间完全独立,不能合并。代码里我会先把所有场景索引对应的变量一次性建好,比如:
soc = sdpvar(n_bus, n_time, n_scenario, 'full'); p_ch = sdpvar(n_bus, n_time, n_scenario, 'full'); p_dis = sdpvar(n_bus, n_time, n_scenario, 'full');这样建变量虽然内存占用大一些,但对求解器来说约束矩阵更规整,不易出错。投资变量是全局的,同一个容量在所有场景下共用。
主问题还有一个细节是投资成本按等年值折算。光伏和风电的寿命一般是20年,储能电池的循环寿命可能只有10年,柴油机更短。代码里用一个折现率r把各设备的一次性投资成本折算成年度等值成本,再和年运行成本相加。这个处理直接影响方案对比的公平性,不能忽视。
4.3 子问题的最坏场景识别
子问题的输入是第一阶段决策x的取值。在实际CCG迭代里,我们不需要求解器显式把x传进去,只需要把主问题求解得到的容量结果赋给子问题里的参数变量,然后重新建模求解即可。
子问题模型包含两部分变量:运行变量y(各时段出力、储能充放电等)和不确定变量u(光伏出力、负荷波动)。目标函数是最大化运行成本,约束包括两步:
- 原内层问题的可行域约束,也就是给定u时y必须满足功率平衡和出力上下限。
- 不确定参数的取值范围约束。
我刚才提到,子问题经过对偶变换后可以把不确定参数按对偶乘子的符号取值。实际操作中代码会写这样一个循环:
for t = 1:n_time if dual_coef(t) >= 0 u_value(t) = u_max(t); else u_value(t) = u_min(t); end end这个dual_coef是从哪里来的?它是对偶变量λ和约束矩阵B相乘后的系数,代表目标函数对不确定参数u(t)的敏感度。如果系数为正,说明这个时段的不确定参数增大时运行成本会上升,那么最坏场景就取它的最大值;反之取最小值。
预算不确定集合的加入会改变这个逻辑,因为这相当于在所有时段里最多允许Γ个参数同时取偏差值,需要额外做一个0/1选择的线性化。代码里建议用YALMIP的binvar定义选择变量,然后加总约束,这样求解器能正确处理。需要注意此时子问题变成混合整数线性规划,单次求解时间会比纯线性规划长一些,但在CCG框架里子问题迭代次数不多,整体还是能接受的。
4.4 CCG主循环的收敛控制
主循环的写法没有太多玄学,但几个细节需要注意。第一是UB和LB的初始值,LB一般可以设一个明确的下界或者直接设0,UB设一个很大的数,然后迭代中逐步收紧。第二是收敛间隙,我习惯设置为0.01或0.005,也就是1%或0.5%的间隙。间隙设太小时迭代次数会明显增加,但对结果的实际改善非常有限,特别是工程项目的容量配置本来就有整定裕度。
第三个细节是前几轮迭代的UB可能波动比较大。这是因为初始场景比较温和,子问题找到的极端场景哪怕还没加入主问题,主问题也能在当前场景集合下找到相对便宜的方案。随着极端场景不断加入,主问题慢慢变紧,方案会趋向稳定。我在代码里会实时绘制LB和UB的收敛曲线,方便观察有没有振荡。
迭代结束后,还要做一次情景回验。这个回验的意思是,把最终得到的容量配置方案放到一个从未参与迭代的随机场景集合里做运行模拟,检查有没有越限。这一步是工程上的双保险,因为在CCG迭代中,你只能保证在“被识别出的”极端场景和它们之间的组合凸包内可行,而真实的连续不确定域是无穷多个场景,CCG的迭代本质是在逐步逼近那个最重要的极端点。理论上C&CG对凸问题能保证有限次收敛到最优,但回验仍然值得做,尤其当你怀疑不确定集合的形状或者线性化近似存在误差时。
5. 算例设计与配置结果分析
5.1 算例参数设置
我拿一个典型的独立微网算例来测这套代码。系统包含光伏电站、风电机组、储能系统、柴油发电机,另外还可以从配电网购电。负荷曲线用夏季典型日数据,峰值负荷约500kW,光伏出力按晴天辐照曲线折算,风速数据取某风场实测日曲线。设备参数大致如下:
| 设备 | 单位投资成本 | 运行维护成本 | 寿命 | 其他说明 |
|---|---|---|---|---|
| 光伏 | 3500元/kW | 0.02元/kWh | 20年 | 弃光率不超过5% |
| 风电 | 6000元/kW | 0.03元/kWh | 20年 | 出力按实际风速曲线 |
| 储能 | 2000元/kWh | 0.05元/kWh循环 | 10年 | 充放电效率95%,SOC范围0.1~0.9 |
| 柴油机 | 1200元/kW | 燃料成本0.6元/kWh | 15年 | 最小出力30%,爬坡约束 |
| 电网交互 | 500元/kW容量费 | 购电价0.8元/kWh | - | 联络线功率上限200kW |
不确定集合的设置,我分别跑了Γ=0、Γ=4和Γ=8三组对比,其中Γ是预算值,代表所有时段里允许同时偏离预测值的时段总数上限。光伏出力预测偏差取±20%,负荷取±10%,风机取±15%。
5.2 迭代收敛情况和配置结果
在默认参数下,CCG算法大约迭代6到9轮收敛到1%间隙以内。前两轮LB和UB差距比较大,第三轮开始差距迅速收窄,最典型的情况是第四轮之后每次UB的下降幅度已经低于总成本的2%。
横向对比三组结果可以看到一个很直观的趋势:
- Γ=0(确定性优化):光伏配了600kW,储能配了800kWh,柴油机只配了200kW,总等年值成本约198万元。
- Γ=4(中等保守):光伏还是600kW左右,储能上升到1000kWh,柴油机维持200kW,电网交互容量从150kW上升到200kW,总成本约215万元。
- Γ=8(高保守):光伏配置增加到680kW,储能直接到1500kWh,柴油机配到300kW,电网容量也顶到上限,总成本约238万元。
有意思的是,光伏容量在Γ增大后并没有简单增加,反而在Γ=4时出现过光伏下降、储能上升的组合。原因不难理解:光伏的出力不确定偏差本身就大,盲目增加光伏容量意味着面对更宽的最坏场景区间,这时候储能反而更划算。所以鲁棒优化做容量配置,最忌讳的就是直觉上“哪个便宜就多配哪个”,必须让算法在不确定性约束下自己权衡。
5.3 典型日运行特征分析
取出Γ=8方案下的最坏场景运行曲线,可以看到几个显著特征:
- 上午9点到下午3点,光伏出力达到预测上界,储能持续充电,同时柴油机以最低出力运行,多余电量全部存入储能。
- 傍晚18点到21点,负荷进入晚高峰,但光伏出力迅速衰减到预测下界,此时储能开始放电,放电功率优先保障负荷,不足部分由柴油机和电网交互补充。
- 夜间22点到凌晨5点,负荷较低,储能不再放电,而是利用电网低谷电价充电,为第二天早高峰做准备。
分段电价下,储能的机会充放电策略体现得很典型。鲁棒优化找到的最坏场景往往不是全天都恶劣,而是集中在光伏出力最大偏差和负荷尖峰叠加的几个时段,算法能精准锁定这些时段说明不确定集合的建模和求解流程是有效的。
6. 调试过程中的常见问题和操作心得
6.1 主问题不可行的典型原因
主问题不可行是CCG迭代中最常见的问题。我排查的时候基本按下面几条主线走:
- 场景集合里包含了一些极端但是物理上无法同时满足的约束组合。比如你给不确定集合设了太大的偏差,同时储能容量又太小,导致所有场景下的SOC递推永远无法满足。这种情况要检查不确定集合的偏差范围和预算值是否合理,可以先用确定性优化跑一版,再把不确定集合逐步放大看临界点在哪。
- 第二阶段约束重复定义了同一组变量。YALMIP不会自动检查变量被重复赋值,如果你在循环里不小心对同一个约束名赋值多次,之前的约束会被静默覆盖,导致模型约束缺失。
- 大M值选得不当。做线性化时M值既要够大,又不能大到让求解器精度爆炸。我一般取相关变量上界的10倍,再检查一下求解器返回的约束违反情况。
6.2 子问题求解缓慢时的加速技巧
子问题如果每次都调用cplex跑MILP,整体时间会非常可观。我的经验是,先检查子问题是否真的是MILP,如果第二阶段全是连续变量,那么对偶变换后得到的max问题可以直接用线性规划求解,速度会快一个数量级。不少同学在建模时习惯性给储能加了一个“是否充放电”的0/1变量,本来是为了避免同时充放电,但其实用约束p_ch * p_dis = 0或者两个方向的功率变量分别设置上限就能实现,整个过程可能完全不需要整数变量。
如果确实需要整数变量,可以用big-M线性化处理双线性项,然后每次迭代时给求解器提供前一轮求解的热启动解。YALMIP的assign和restore功能可以用来保存上次的变量值,cplex自身也支持warm start。实测下来热启动能节省大约30%到40%的子问题求解时间。
6.3 参数灵敏度分析的快速实现方法
调参在鲁棒优化里是家常便饭。我的建议是把不确定集合参数、设备成本参数、负荷曲线全部封装成结构体,作为主函数的输入。这样你在脚本里改一个字段就能重新跑一整轮分析,特别适合批量做灵敏度扫描。
实际测试结果里,储能成本和光伏成本对配置结果最敏感。如果储能成本从2000元/kWh下降到1500元/kWh,Γ=8方案里的储能配置会从1500kWh涨到接近2300kWh,光伏反而下降100kW。这说明系统在鲁棒性需求下更依赖储能来吸收不确定性,而不是单纯扩大新能源容量。这也是为什么高比例新能源微网里储能几乎成了标配。
7. 边界条件和扩展方向
两阶段鲁棒优化这套框架并不局限于微网容量配置。同一个代码骨架,改一下不确定集合的定义和第二阶段的操作约束,就能用来做配电网网架规划、综合能源系统容量配置、备用容量优化等。目前已经有不少研究把CCG方法推广到多能互补系统,比如加入热、气等多种能源形式,虽然第二阶段运行调度的复杂度会成倍上升,但基本思想是一致的。
还有两个值得关注的扩展方向:一是和数据驱动相结合,不再用简单的盒式集合,而是从历史数据中提取场景并建立数据驱动的模糊集,用分布鲁棒优化替代传统鲁棒优化,这样能在不牺牲太多保守度的情况下提高经济性;二是和强化学习结合,用鲁棒优化生成容量方案,用强化学习做实时运行策略,两阶段问题的第二层运行决策可选择的求解方法会更灵活。
就我的使用体验来说,MATLAB这套CCG代码最适合的还是作为研究原型验证工具。做工程项目的同学,后续可以考虑用Julia的JuMP或者Python的Pyomo把同一套模型复现一遍,因为这两个生态在处理大规模不确定优化时性能更好,也不存在和操作系统绑定的问题。不过如果只是想快速验证算法思路、跑通算例、出曲线交差或者发论文,MATLAB足够顺手了。
最后再分享一个代码实践的小经验:不要把所有内容都堆在一个脚本文件里。按上面的模块划分,数据、建模、求解、后处理各放各的文件,接口用结构体传参,后续调参会轻松很多。CCG这类型迭代算法,代码的中间输出特别重要,每轮迭代的结果、间隙、场景信息建议全部保存到workspace变量或者输出到Excel,方便回溯。写代码的时候多留一步日志,排查问题的时候能省出几倍的时间。