两阶段鲁棒优化在微网容量配置里这几年是真的火，但也是真的让人头大。我在做微网规划项目的时候，一开始用确定性方法算出来的容量方案，拿到实际运行数据里一测，不是光伏出力波动稍微大一点就出现切负荷，就是柴油机组被频繁启停折腾得寿命骤降。后来把两阶段鲁棒优化框架搭起来，才算是把"容量决策"和"运行调度"这两层问题真正捏合到了一起。

这篇就围绕我的Matlab实现过程，把从建模、求解到代码落地的一系列关键细节讲清楚。项目主场景是含光伏、风电、柴油发电机和储能电池的并网型微网，目标是给出既能应对最恶劣风光出力场景、又不会过度浪费投资的电源容量配置方案。如果你想复现一个可运行的鲁棒容量配置Demo，或者正在纠结怎么把C&CG算法落实到代码里，这篇应该能省下你不少调研时间。

1. 微网容量配置为什么要引入两阶段鲁棒优化

很多刚接触这个领域的人会有个疑问：传统的确定性容量配置算得好好的，为什么非要用鲁棒优化？我在项目里体会最深的一点是，微网容量配置本质上是个"决策在前、验证在后"的问题。你要在规划阶段把光伏板、风机、柴油机、储能的容量定下来，但真正运行的时候，光照、风速、负荷都是不确定的。如果只用某几个典型场景的期望值来做决策，一旦实际风光出力落在比较差的区间，容量配置就可能出现系统性缺电。

1.1 风光出力的不确定性与容量配置的本质困境

微网里的光伏出力和风电出力从来都不是恒定值。某天早上太阳好，光伏能做到额定容量的90%以上，到了下午云层一遮，瞬时出力能掉到20%。风速更是毫无规律的，我给项目做数据统计的时候，某地冬季风电出力的标准差能到额定容量的30%以上。这种情况下，如果你在规划模型里把风光出力当成确定参数，算出来的容量配置一定是"顾头不顾尾"的。

容量配置还有一个麻烦是"双层时间尺度耦合"。投资决策是年尺度甚至二十年尺度的，运行调度却是分钟级到小时级的。你不能只算设备装多少，还得验证装完之后能不能在每个典型时刻满足功率平衡和备用要求。这就意味着容量配置本质上是一个带运行约束的优化问题，计算复杂度比单纯的设备选型高很多。确定性方法之所以不够用，是因为它把不确定参数全部用期望值或者典型场景替代了，本质上是拿一个"平均情况"去代表所有可能出现的运行场景，自然扛不住恶劣场景的冲击。

1.2 随机优化与鲁棒优化：一个硬币的两面

处理不确定性的主流方法有两类。随机规划需要假设不确定参数的概率分布已知，比如光伏出力的预测误差服从正态分布，然后用场景树去做期望值优化。它的好处是结果经济性通常比较好，不会为极端场景付出过高的容量成本，但坏处是你必须对概率分布有足够置信度，而且场景数量一多，求解规模会爆炸。

鲁棒优化走的是另一条路。我不需要知道光伏出力波动的精确概率分布，只需要定义它可能波动的范围，然后在一个不确定集合里面搜索最坏情况。最坏情况下系统也能安全运行，那么其他场景大概率也能扛住。容量配置这种决策问题，设备的安装容量一旦买入就是沉没成本，你想反悔都不行，所以用保守的鲁棒策略反而有天然的适配性。

我在项目里对比过随机规划和鲁棒优化的结果，随机规划给出的光伏容量会比鲁棒方案低10%左右，总成本也低一些，但一旦把实际最差出力场景套进去，随机规划方案会出现最大20%的甩负荷量。对微网来说，供电可靠性往往比一两年的成本差额重要得多，鲁棒优化的价值就在这里。

1.3 两阶段决策结构：容量先行，调度随动

两阶段鲁棒优化的"两阶段"不是指时间上的两个步骤，而是指决策变量的两个层次。第一阶段变量是"现在就要定下来"的决策，对应到容量配置里就是光伏装机容量、风机装机容量、储能额定容量和功率、柴油机台数。这些变量一旦确定，在后续所有场景里都不能变。第二阶段变量是"看到不确定性实现之后才能做"的决策，对应的是每个运行时刻的机组出力、储能充放电功率、弃风弃光量、切负荷量等运行调度量。

用数学语言来说就是：

min（投资成本 + 最坏情况下的运行成本）

第一阶段决定容量，然后自然就产生一个问题：容量定了之后，在最坏的风光出力场景下，运行成本是多少、会不会出现失负荷。这就是内层的max-min问题。外层求容量最优，内层搜最坏场景并给出相应调度成本，两者嵌套在一起，构成了一个典型的两阶段鲁棒优化模型。

我在写代码的时候，一开始试图把两个阶段合并成一个单层混合整数规划来求解，结果变量数直接爆炸，根本跑不动。后来老老实实换成C&CG迭代求解框架才算是找到了正路。

2. 数学模型搭建：从投资层到运行层的双层联动

模型是整篇代码的地基，建模建歪了后面所有求解器调试都是在浪费时间。我这边以典型的光-风-柴-储微网为例，把核心变量、目标函数、约束条件和不确定集合逐个拆开讲。你可以直接拿这套模型作为自己项目的骨架，根据实际情况增加或删减设备类型。

2.1 目标函数：投资成本与运行成本的"双时间尺度"耦合

两阶段鲁棒优化的目标函数分两层。第一阶段目标函数是设备投资成本的年化值，包含光伏、风机、柴油机组、储能系统四项。每个设备都要考虑单位容量投资成本、等年值系数，以及设备的使用寿命。等年值系数就是把一次性投资按照折现率和使用年限换算成每年摊销金额，公式是：

r × (1+r)^n / ((1+r)^n - 1)

其中r是折现率，n是设备寿命。我用的是8%折现率，光伏和风机按20年寿命算，储能按10年，柴油机按15年。这套参数也是文献里常用的，算出来结果比较有参考价值。

第二阶段目标函数是运行成本，在给定容量方案和给定风光出力场景下，包含柴油机燃料成本、启停成本、储能老化成本、弃风弃光惩罚成本，还有切负荷惩罚成本。切负荷惩罚要设得足够高，我一般设成10元/kWh以上，否则优化器会倾向于用切负荷代替昂贵的燃料成本或储能耗损。

两层目标函数通过一个"最坏场景"连接起来。第二阶段运行成本不是某个固定值，而是在不确定性集合内所有可能场景中的最大值。整体形式如下：

min（投资年化成本 + max（运行成本））

这个max内层是模型里最核心的部分，也是求解难度最大的来源。

2.2 约束条件：功率平衡、储能电量和机组出力

约束条件同样分两个阶段。第一阶段的约束主要是设备容量变量的非负约束和上限约束，比如安装面积决定了光伏容量上限，变压器容量决定了总装机上限。这些比较简单，不多说。

第二阶段的运行约束就细了。首先每个时刻都要满足功率平衡：

光伏出力 + 风电出力 + 柴油机出力 + 储能放电功率 + 外网购电 = 负荷 + 储能充电功率 + 弃风弃光量

注意这里我特意加上了弃风弃光量和切负荷量这两个松弛变量。没有它们模型容易无解，因为当风光出力特别高而负荷很低时，系统没法把多余的电量完全吸收，必须允许弃掉一部分；反过来，当风光出力特别低时，如果所有机组都开满还是不够，就需要切负荷来保证功率平衡。松弛变量进入目标函数时带惩罚系数，这样优化器不会滥用它们。

储能约束是整个模型里最需要小心的地方。储能SOC（荷电状态）的时序递推关系是：

SOC(t+1) = SOC(t) + η_c × P_ch(t) - P_dis(t) / η_d

其中η_c和η_d分别是充放电效率。充放电功率不能同时为正，这个我会通过引入二进制变量来控制，也可以直接把充放电功率设成双向变量，然后用它的正负和效率函数来处理。

柴油机组需要建最小启停时间约束吗？理论上需要，但实际建模时如果调度周期是24小时，引入最小启停时间会带来大量二进制变量，大大增加求解难度。我在项目里折中了一下，只加了机组出力上下限约束和爬坡约束，没有加最小启停时间约束，这样既保留了运行调度的关键物理特性，又不至于让模型彻底丧失可求解性。

2.3 不确定集合的建模：盒式、椭球和预算约束

两阶段鲁棒优化的核心在于不确定集合的构建。最常见的集合形式是盒式集合，就是每个不确定参数各自在某个区间内独立波动。比如光伏出力P_pv可以在预测值上下偏差±ΔP_pv_max的范围内变动，风电也是类似。

盒式集合最简单，但有一个问题：它允许所有不确定参数同时达到最坏取值。现实中不会出现所有光伏组件同时出力最低又所有风机同时停转的情况，更合理的假设是最多只有一部分不确定参数同时偏离到极端值。这就需要引入预算约束：

Σ |ΔP_i| / ΔP_i_max ≤ Γ

其中Γ是鲁棒控制参数，它限制了所有不确定参数总偏离程度的上限。Γ取0就是完全确定性问题，Γ取最大允许值就是盒式集合的最保守情况，Γ取中间值则可以在保守性和经济性之间调节。

我实际用的就是这个带预算约束的盒式集合，具体是在每个时段内限制风光出力偏差的累计程度。也有文献用椭球集合，它更适合不确定参数之间存在相关性的场景，但椭球集合的求解需要二阶锥约束，在某些求解器组合下处理起来比较麻烦。做容量配置这种偏工程的项目，带预算约束的盒式集合基本够用了，而且对Matlab代码实现也更友好。

2.4 为什么预算约束能控制保守程度

很多刚上手的人对Γ的理解比较抽象。我用一个生活化的类比说明：假设你要出门旅行，正常情况下路上会有交通拥堵、天气延误、出租车排队这些小意外。如果假设所有最坏情况同时发生，你会准备一条几乎不可能用完的巨额应急资金，这就是盒式集合的过度保守。但如果假设最多只允许两个坏情况同时发生，你准备的应急资金就会合理得多，这就是预算约束的价值。

在微网容量配置里，Γ控制的正是"最坏场景有多坏"。Γ越小，你准备应对的最坏场景越温和，算出来的容量成本越低，但实际碰到更恶劣场景时的失负荷概率也越高。我在项目里会做一组Γ扫描实验，比如从1到最大偏差量，观察容量配置结果和总成本的变化曲线，然后选择一个成本增长斜率突变的点作为推荐值。

我自己的经验是，Γ取预算最大值的40%到60%通常是比较合理的区间，这个范围内能在可靠性保证和投资成本之间找到不错的性价比。当然不同地区的光照和风速相关性不一样，建议你拿到数据后先做统计，看看本地区风光出力同时波动的实际规律再定。

3. 求解框架：列与约束生成算法（C&CG）的原理与实现

模型建完之后，真正的技术难点在求解。两阶段鲁棒优化是个min-max结构，不能直接用现成的优化求解器求解，必须用迭代算法。我在项目里用的是列与约束生成算法，也就是C&CG，它比传统的Benders分解在这个问题上收敛快得多，代码实现也相对清晰。这一节是我自己从零搭完整套求解流程之后总结出来的，能帮你少走不少弯路。

3.1 主问题MP与子问题SP的博弈结构

C&CG的核心思想很简单：把原问题拆成主问题和子问题，两者互相博弈，交替求解直到收敛。

主问题MP是一个有限场景的优化问题。原始的不确定集合包含无穷多个可能的场景，不可能全部列出来，但主问题里只保留已经被识别出的有限个关键场景，求一个在这些场景下都可行且成本最低的容量配置和调度方案。主问题求解完，得到容量方案x和总成本下界LB。

子问题SP的输入就是主问题算出来的容量方案x。子问题要在整个不确定集合内搜索一个最坏场景，让给定x下的运行成本最高，同时给出这个最坏场景下的最优调度方案。子问题求解完，得到最坏场景和对应的运行成本，把它和当前容量方案的投资成本相加，得到总成本上界UB。

如果UB和LB的差距缩小到预设的容差内，说明当前容量方案在最坏场景下的成本和主问题内场景的成本已经基本一致，迭代就可以终止了。如果差距还大，就把子问题找到的最坏场景作为新场景加到主问题里，重新求解，如此反复。

这个过程有点像一个谈判：主问题先出一个容量方案，子问题立刻指出"你这个方案在我手上这个最坏场景里成本很高"，主问题收到反馈后修正方案，子问题再找更坏的场景，直到双方达成一致。

3.2 子问题的最坏场景搜索：对偶、KKT与Big-M

子问题的结构是max-min，外层要在不确定集合内找最坏的风光出力场景，内层要在最坏场景下求最小运行成本。这不是一个可以直接求解的标准数学规划问题，需要处理。

最常用的处理手法是把内层min问题通过拉格朗日对偶转成max形式，然后和外层max合并，变成一个单层max问题。这个单层问题里目标函数和约束条件中会出现不确定变量u和调度变量y的乘积项。按照标准鲁棒对偶理论，这类乘积项可以通过引入对偶变量来处理，最终得到的是一个可以在求解器里求解的优化问题。

但这里有个坎：如果内层问题存在整数变量，比如储能充放电状态的二进制变量，那么对偶理论就不直接适用了。不是不能用KKT条件或Big-M线性化来处理，但处理难度会上升。我的做法是尽量让内层问题保持为线性规划可解的形式，具体办法是：储能充放电状态不用二进制变量区分，而是把充放电功率拆成两个非负变量，并用惩罚项和效率约束来隐式约束不同时充放。这样内层问题保持线性，对偶转换就干净了。

如果必须保留整数变量怎么办？可以把内层问题用KKT条件替换，再通过Big-M把互补松弛条件线性化，得到一个混合整数规划子问题。这种方法可行，但Big-M参数需要仔细选择，否则数值稳定性会很差。我在测试中发现，M值设成负荷最大值的三到五倍通常能兼顾收敛性和数值稳定。

3.3 迭代收敛判定

C&CG的收敛判定很简单，就是持续比较上下界。

初始状态下，LB可以设成负无穷，UB设成正无穷。每轮迭代，主问题求得LB，子问题求得UB。当(UB-LB)/UB小于预设阈值，比如0.01或者0.001，就认为收敛。阈值设太大会导致容量方案偏离真实最坏场景太远，设太小则迭代次数暴增，求解时间从几分钟变成几小时都是可能的。

我在代码里用的是相对间隙1e-3，配合一个迭代次数上限30次。实测下来，对带预测不确定性的容量配置问题，一般15到25轮之间能收敛。如果超过30次还没收敛，说明模型某个地方大概率有数值问题，而不是单纯收敛慢。这时候最有效的排查方式是检查子问题对偶变量的尺度是否相差过大，以及Big-M参数是否设得离谱。

另外有个小技巧，主问题每一轮迭代都会积累新的场景和变量，场景变量越来越多，主问题规模越来越大。所以必须在主问题里把每一轮迭代生成的最坏场景都保留下来，不能只保留最新一轮的。这不是写代码时容易漏掉的地方，而且一开始使用MPCC工具时特别容易误以为保留最新场景就够了。

3.4 C&CG与Benders分解的取舍

Benders分解和C&CG都能求解两阶段鲁棒问题，但我在实际对比中明显感觉到C&CG在容量配置这个场景下更优。原因在于C&CG每次迭代向主问题添加的是完整的场景变量和约束，主问题的约束越来越贴近真实的可行域，收敛速度显著更快；Benders分解添加的是割平面近似，收敛比较慢，尤其当最优解落在可行域内部的时候会震荡。

Benders分解的一个优势是主问题规模增长速度慢，内存占用较小。但对于微网容量配置这种中小规模优化问题，C&CG增加的那点规模根本无所谓。我做过一个含24时段、5种设备类型的测试算例，C&CG大概25轮收敛，Benders跑了60轮还在震荡，结论非常明显。

如果你用的是成熟的鲁棒优化工具包，可能不需要自己实现C&CG。但我还是强烈建议你至少手动实现一遍C&CG流程，因为只有拆开理解主问题和子问题的博弈关系，才能在处理模型的特殊约束时知道该改哪里。

4. Matlab实现的关键代码逻辑与工程细节

Matlab生态下的组合是YALMIP+Gurobi或者YALMIP+Cplex。把C&CG算法用YALMIP表达出来，好处是可以用数学表达式直接写约束，可读性高，调试时容易定位问题。坏处是YALMIP本身是建模层的解释器，如果模型规模大，它的符号运算开销会拖慢整体效率。好在容量配置问题规模不大，YALMIP完全够用。

4.1 工具链选型：YALMIP+求解器的搭配逻辑

YALMIP负责把优化模型翻译成求解器能理解的标准形式。线性规划、二次规划、混合整数规划都有对应的求解器支持。在我的项目里，主问题是一个混合整数线性规划，因为储能充放电功率拆成两个变量后不再需要二进制变量，但柴油机启停和容量离散选择会产生0-1变量。子问题经过对偶转换后是一个线性规划或者带二阶锥约束的问题。

选择具体求解器时，Gurobi和Cplex的线性规划、混合整数规划性能都非常强，在中小规模问题上差距不大。如果学生党没有商业求解器授权，用开源的Cbc、GLPK也能跑通流程，只是速度慢一些。我在机器上装的是Gurobi，授权的事就不多提了，你自己根据条件选择就行。

另外需要提醒的是，YALMIP版本和求解器版本的兼容性经常出问题。我遇到过YALMIP调用求解器返回空解的情况，最后排查下来是YALMIP版本太旧，不识别求解器的新参数。建议用较新的YALMIP版本，并且求解器也保持更新。

4.2 代码骨架与数据结构

整个项目的代码结构我分成三块：参数初始化脚本、C&CG主循环函数、结果后处理脚本。

参数初始化脚本里定义所有系统参数，包括各类设备的成本参数、效率参数、预测出力曲线、负荷曲线、不确定集合的偏差范围和预算Γ、折现率、寿命年限等。一些重要的心得是：把参数集中放到一个结构体里，以sysParam命名，后面每个函数都传入这个结构体，比全局变量好用得多。

数据这块最费时间的是准备预测出力曲线。我用的方法是，把历史光伏出力和风电出力按小时取平均，得到基准预测曲线，然后根据历史数据的标准差确定每个时段的不确定偏差范围。如果你手头没有历史数据，可以用合成方法生成，比如用正弦曲线拟合日照规律，再叠加上下偏差，但这样出来的结果会偏理想，不一定能反映真实波动特征。

C&CG主循环函数是核心。函数先初始化UB和LB，然后循环调用主问题构造函数和子问题构造函数。主问题构造函数接收已经发现的所有最坏场景集合，用YALMIP定义容量变量和对应各场景的运行变量，求解后返回容量值、目标函数值和割平面信息。子问题构造函数接收主问题给的容量值，构建对偶转换后的最坏场景搜索模型，求解后返回新的最坏场景和对应运行成本。

结果后处理脚本把每个设备的最优容量、总成本、各场景的调度方案等结果输出成表格和曲线，方便比对。

4.3 核心代码段演示

我先展示主问题的核心代码逻辑，省略设备参数定义部分，只保留建模核心部分：

function [x_opt, obj_mp, mp_model] = build_MP(scenarios, sysParam) % 定义第一阶段容量变量 x_pv = sdpvar(1, 1, 'full'); x_wt = sdpvar(1, 1, 'full'); x_ess = sdpvar(1, 1, 'full'); x_dg = sdpvar(1, 1, 'full'); % 投资成本年化值 inv_cost = sysParam.cost_pv * x_pv * sysParam.crf_pv + ... sysParam.cost_wt * x_wt * sysParam.crf_wt + ... sysParam.cost_ess * x_ess * sysParam.crf_ess + ... sysParam.cost_dg * x_dg * sysParam.crf_dg; % 循环每个场景定义运行变量和约束 run_cost = 0; constraints = []; for k = 1:length(scenarios) % 场景k的风光出力和负荷数据 p_pv_k = scenarios(k).p_pv; p_wt_k = scenarios(k).p_wt; p_load_k = scenarios(k).p_load; % 第二阶段运行变量 p_dg = sdpvar(24, 1, 'full'); p_ess_ch = sdpvar(24, 1, 'full'); p_ess_dis = sdpvar(24, 1, 'full'); p_curtail = sdpvar(24, 1, 'full'); p_shed = sdpvar(24, 1, 'full'); soc = sdpvar(25, 1, 'full'); u_dg = binvar(24, 1); % 柴油机启停状态 constraints = [constraints, ... 0 <= p_dg <= sysParam.p_dg_max * x_dg, ... 0 <= p_ess_ch <= sysParam.p_ess_rated * x_ess, ... 0 <= p_ess_dis <= sysParam.p_ess_rated * x_ess, ... 0 <= p_curtail <= p_pv_k + p_wt_k, ... 0 <= p_shed <= p_load_k, ... soc(1) == sysParam.soc_init, ... soc(25) == sysParam.soc_init, ... soc(2:25) == soc(1:24) + sysParam.eta_ch * p_ess_ch - p_ess_dis / sysParam.eta_dis, ... sysParam.soc_min <= soc(2:25) <= sysParam.soc_max, ... p_pv_k + p_wt_k + p_dg + p_ess_dis + p_shed + p_curtail == p_load_k + p_ess_ch]; % 运行成本 run_cost = run_cost + sum(sysParam.fuel_cost * p_dg) + ... sum(sysParam.start_cost * max(0, diff([0; u_dg]))) + ... sysParam.penalty_shed * sum(p_shed) + ... sysParam.penalty_curtail * sum(p_curtail) + ... sysParam.ess_degrad_cost * sum(p_ess_ch + p_ess_dis); end % 求解主问题 optimize(constraints, inv_cost + run_cost, sdpsettings('solver', 'gurobi', 'verbose', 0)); x_opt = value([x_pv, x_wt, x_ess, x_dg]); obj_mp = value(inv_cost + run_cost); end

注意这里max(0, diff(...))这种表达在YALMIP里不是标准线性形式，实际项目里我会引入辅助变量来表示柴油机启停成本，这里写成简单形式是为了展示逻辑。实际代码里应该用线性化约束：

z_t >= u_t - u_{t-1}; z_t >= 0

来表示启停动作。

4.4 子问题的对偶转换代码实现

子问题实现比主问题麻烦，需要对内层调度问题做对偶转换。但YALMIP有一个可以偷懒的路径：用YALMIP的robust框架或者uncertain声明来直接建模鲁棒优化。我在测试初期试过这个功能，不过在带预算约束的不确定集合下，YALMIP自带的鲁棒建模功能生成的对偶形式经常不是我想要的，而且耦合了自己写的双层结构时容易报错。

所以最终我选择手动写对偶转换。做法是先定义内层调度问题的原始变量和对偶变量，然后用符号推导写出对偶问题的约束和目标函数，再与外层max合并。这个过程最容易出错，我建议你在测试小规模问题时先用一个3时段的简化模型验证对偶是否正确，再扩展到24时段。

对偶转换完成后，子问题的目标函数会呈现类似下面的结构：

function [u_worst, obj_sp, sp_model] = build_SP(x_opt, sysParam) % 不确定变量 delta_pv = sdpvar(24, 1, 'full'); delta_wt = sdpvar(24, 1, 'full'); % 不确定范围约束 constraints = [... -sysParam.delta_pv_max <= delta_pv <= sysParam.delta_pv_max, ... -sysParam.delta_wt_max <= delta_wt <= sysParam.delta_wt_max]; % 预算约束 constraints = [constraints, ... sum(abs(delta_pv) ./ sysParam.delta_pv_max) + ... sum(abs(delta_wt) ./ sysParam.delta_wt_max) <= sysParam.Gamma]; % 对偶变量和线性化后子问题目标 % 具体形式取决于原始内层约束的系数矩阵，需要根据你的模型推导 % 这里用简化示意 obj_sp = ...; % 对偶目标函数，包含x_opt和delta的耦合项 optimize(constraints, -obj_sp, sdpsettings('solver', 'gurobi', 'verbose', 0)); u_worst.p_pv = sysParam.p_pv_base + value(delta_pv); u_worst.p_wt = sysParam.p_wt_base + value(delta_wt); obj_sp = value(obj_sp); end

绝对值项在YALMIP里处理需要额外变量，实际项目中我用的是引入正负偏差变量的方式，把delta_pv拆成delta_pv_pos和delta_pv_neg，然后用它们的和来构建预算约束。这样整体保持线性约束的表达形式。

4.5 求解器参数和数值预处理

YALMIP求解时的sdpsettings参数里有几个值得调。一个是数值容差，Gurobi默认的绝对间隙和相对间隙设置比较保守，我一般会把相对间隙设成1e-4或者更高，避免主问题和子问题在边界上反复震荡。另一个是线程数，如果模型规模不大，设成单线程反而更稳定，多线程在小模型上的收益不明显，还容易引入结果不确定性。

数值预处理方面，建议把所有成本系数和功率量纲统一。我习惯把所有功率用kW，电量用kWh，成本用元，这样数量级差异不会太大。如果有些系数相差超过1e4倍，求解器的数值稳定性会明显下降，出现莫名其妙的不可行或者无界警告。检查数值尺度是否合理，是调试C&CG循环的第一步。

5. 实测中的坑与调参经验

这部分是最干货的部分，都是我在代码跑通和项目复现过程中真正踩过的坑。如果你是自己动手写代码，大概率会撞上其中几个。

5.1 双线性项的处理：别偷懒

C&CG子问题里，不确定变量u和对偶变量λ会以乘积形式出现在目标函数中，这是max-min结构带来的天然产物。YALMIP可以直接建立两个sdpvar变量的乘积，但这样得到的模型是双线性的，Gurobi根本解不了。

我在处理这个双线性项时的经验是：如果是连续变量和连续变量的乘积，根据不确定变量的边界值，可以通过大M法或者McCormick包络线性近似，或者直接利用盒子约束的对偶性质消去变量；如果是二进制变量和连续变量的乘积，则要先引入中间变量，再用大M法线性化。这里有个容易踩的坑：McCormick包络是松弛近似，得到的最坏场景不一定是最优最坏场景，保守性会受影响。所以能通过严格对偶变换直接消去双线性项的情况下，尽量用严格方法。

5.2 迭代不收敛的排查路线

如果C&CG迭代超过40轮还在上下界之间震荡，不要盲目加大迭代上限，这是一条错误的排查路径。更有效的做法是检查以下几点：

第一，检查子问题是否真的求到了最坏场景。把子问题算出来的风光场景和目标函数值记录下来，人工验证一下在这个场景下，主问题给出的容量方案是否真的会产生这么高的运行成本。如果手工验算成本比子问题给出的低，说明子问题的对偶转换或者双线性项处理错了。

第二，检查主问题是否保留了所有历史场景。我在初版代码里犯过这个错，主问题只保留最新场景，结果每次迭代都在不同的场景之间跳来跳去，上下界永远收敛不到一起。把场景保留改成累积式之后，大约第八轮就开始稳步收敛了。

第三，检查数值容差。上下界差距缩小到一定程度就不动了，可能是求解器精度限制。试着把求解器的相对间隙收紧一点，或者把目标函数系数统一缩放1000倍，排除数值尺度问题。

5.3 保守性-经济性平衡的调参

Γ的值直接决定了容量配置结果的保守程度。我在项目里做了个实验，把Γ从0逐步加到最大预算值，观察光伏、风机、储能容量的变化。

在没有预算约束时，即Γ为0，系统完全不考虑不确定性，光伏容量大约是600kW，储能配置也会比较小，总成本最低。随着Γ逐步增大，系统为了应对更恶劣的风光场景，光伏和风机容量反而可能减小，因为不可靠的电源容量在极端场景下形同虚设，柴油机和储能的容量则稳步上升。在Γ达到最大值时，储能容量可能比Γ为0时高出50%以上，柴油机容量也显著增加。

这个趋势说明，如果某个地区的风光资源波动特别大，单纯加大光伏和风机容量并不能提升供电可靠性，反而应该配置灵活调节资源，如储能和柴油机。在做项目汇报的时候，我会把这张Γ-容量曲线直接放进去，决策者看到之后更容易理解保守程度对投资的影响。

5.4 24时段扩展到8760小时的加速思路

容量配置严格来说应该用全年8760小时的时序数据做校核，但两阶段鲁棒优化的子问题本身是24时段，直接扩展到8760会让变量规模爆炸。我试过直接扩，结果主问题单轮求解就要半小时以上，迭代根本走不动。

工程上我最后采用的方法是典型日法加季节性分解。把一年分成四个季节，每个季节选一个代表日，然后对四个典型日分别做容量配置，取各设备容量的包络值作为最终结果。这样虽然会牺牲一点精度，但计算时间能控制在十几分钟内，而且结果的可解释性更强。

更进阶的做法是用Benders分解和C&CG的嵌套结构，把8760小时的运行调度问题拆到多个时段块分别求解。这个方向我在后续迭代里验证过，速度和精度都比典型日法好，但代码复杂度明显上升。如果只是做一个初步的可行性研究，典型日法已经够用了。

5.5 关于Matlab代码复用和扩展的几个建议

代码结构上，强烈建议把C&CG主循环、主问题构造、子问题构造、参数设置、结果输出拆成独立函数，每个函数用清晰注释说明输入输出。我一开始图省事写了个巨型脚本，调试时改一个参数要拖半天，后来重构成函数化结构后，加新设备、改不确定集合都很方便。

不确定参数的波动范围设置要基于实际数据。如果项目手头暂时没有历史数据，至少也要参考行业报告或仿真数据，不要随意设一个对称区间，那样算出来的配置方案缺乏实际参考价值。

如果你想在这个基础上加新的不确定性来源，比如负荷不确定性或者电价不确定性，只需要在不确定变量列表里增加对应变量，在预算约束里加入对应的偏离项即可，整体框架不用大改。C&CG算法的兼容性很好，这也是它在工程界受欢迎的原因之一。

最后提一个我踩过最多次的坑：每次迭代开始前，记得清空YALMIP的变量缓存和求解器的工作区，否则旧模型的变量约束会残留，造成模型规模膨胀和求解偏差。YALMIP里sdpvar对象重复定义不会直接报错，你以为是在更新模型，实际上是在叠加约束，这个问题排查起来相当隐蔽。