1. 项目缘起这不是“套模型”是一场和不确定性较劲的实战做微电网优化调度的人应该都有同感你辛辛苦苦搭了一个确定性经济调度模型把风光预测值往里面一塞算出来一套调度计划以为万事大吉。结果第二天光伏被云遮了半小时风机出力跟预测差了快两成实时平衡市场电价一波动原本算好的“最优”方案瞬间变成了一笔亏钱方案。这种挫败感做过工程的人都能懂。我这次复现的两阶段鲁棒优化经济调度方法解决的就是这个问题。它的核心思路是不在预测值上赌运气而是假设不确定性在某个给定范围内任意变化找一套在所有可能场景下都能保证可行、并且最坏情况下总成本最优的调度方案。说白了就是用一点点经济性代价换取系统运行的安全性边界。这套方法在当前微电网、主动配电网、综合能源系统的研究中是妥妥的热点方向也是工业界从“计划调度”走向“韧性调度”绕不开的一步。这篇博文不是给论文做摘要而是把从数学模型、算法设计、代码复现到结果验证的完整链路记录下来。适合三类人看一是正在复现相关论文、卡在两阶段鲁棒建模细节里的研究生二是想在企业微电网尤其是矿山、工业园区这类对供电可靠性要求苛刻的场景落地鲁棒调度的工程师三是刚开始接触鲁棒优化、想从零理解“列与约束生成”到底在干什么的入门者。如果你是第三类建议把第一节和第二节反复看两遍后面再跟着代码实操。2. 两阶段鲁棒优化到底在干什么先把这个模型掰开了理解2.1 两阶段的结构先定“现在”再应对“未来”两阶段鲁棒优化的名字听起来高大上但拆开看就是一个“现在做决定、未来做调整”的决策框架。拿微电网来举例调度员在日前比如上午10点就需要确定次日的机组启停、与主网的购售电计划这是第一阶段决策特点是决策时不知道确切的分布式电源出力和负荷。等第二天实际运行到某个时刻风光出力、负荷的真实数据暴露了一部分调度员需要在第一阶段决策的基础上对机组出力爬坡、电池充放电功率、切负荷量做实时调整来保证系统功率平衡这是第二阶段决策。这里有几个容易绕晕的点我展开说一下第一阶段决策也叫“here-and-now”决策特点是决策变量在不确定性实现前就必须敲定具有不可更改性。比如火电机组的启停状态你不能说“等风来了再决定启不启机”。第二阶段决策叫“wait-and-see”决策也就是看到不确定性实现值之后再做的适应性调整。它存在的前提是系统有灵活性资源比如储能、可调负荷、机组爬坡能力。鲁棒优化和随机规划最本质的区别就在这里随机规划给不确定性变量指定概率分布然后用期望值目标鲁棒优化不给分布只给定一个不确定集合然后找最坏情况下的最优解。如果你从来没接触过这类模型可以拿买房做类比。第一阶段决策是“签贷款合同”这时你锁定了贷款方案第二阶段决策是“交月供”如果利率上浮不确定性实现你需要调整每月还款结构来保证不违约。鲁棒优化的思路就是提前假设利率可能在某个区间内波动然后找一个无论利率怎么变都不会违约、并且最坏情况下总还款额最低的贷款方案。2.2 模型表达式用数学语言把这个调度问题锁死我们今天复现的微电网拓扑结构是典型的交流微电网带储能系统包含燃气轮机、储能电池、风机、光伏、固定负荷和可调负荷同时和配电网有联络线可以双向购售电。这里我把核心数学模型写出来注意这里用的不是论文里那种缩写极多的套路而是可以直接翻译成代码的写法。目标函数是所有阶段的总成本最小包括第一阶段成本燃气轮机启停成本、日前与主网的购电成本。第二阶段成本燃气轮机燃料成本、某些可调负荷的响应补偿成本、弃风弃光惩罚成本、最坏场景下的切负荷惩罚成本。写成标准的两阶段鲁棒优化形式就是min c1^T y max u∈U min x(z,u) c2^T xs.t. Ay ≥ b 这是第一阶段约束比如机组启停逻辑、最小启停时间约束Bx ≥ d - Cu - Ey 这是第二阶段约束比如功率平衡、机组出力上下限、储能SOC递推公式u ∈ U U就是我们设定的不确定集合。这里必须强调一个理解上的关键卡点整个模型是一个“min-max-min”三层结构。内层的min是第二阶段优化问题max是找到一个最坏的不确定性实现u让总成本最大外层的min是选择第一阶段决策变量y让最坏情况下的成本最小。这个结构直接求解是不可能的因为内层的max会把可行域搞得不可凸所以必须用算法去迭代逼近。目前主流的解法是两个列与约束生成算法也是我们这次复现用的方法。它的思想是先固定一个不确定场景求解主问题得到一个下界然后把主问题的解代入子问题求解一个max-min问题找到真正最坏的不确定性场景再把这个场景作为新的约束列加入到主问题里反复迭代直到上下界收敛。Benders对偶分解思想是把max-min子问题通过对偶变换转化为单层max问题然后用割平面去逼近。这个方法的缺点是当第二阶段是整数变量时不好用。选择CCG而不是Benders原因很直接CCG在处理含整数变量比如可调负荷是否参与响应的0-1状态时不需要做主问题的多次可行割回代收敛速度要快得多。在实际的微电网调度中储能充放电状态、可调负荷的投入状态这类整数变量几乎必然出现所以CCG是更稳妥的选择。2.3 不确定集合的设计决定鲁棒性和经济性的天平很多人在复现论文时把精力全花在求解算法上却忽视了在最前面就要做的一个决定——不确定集合长什么样。不确定集合就是给风光出力和负荷预设一个波动区间形式直接决定了模型的保守程度。常见的有三种盒式集合。最简单每个不确定性变量独立在一个区间内变化。缺点是它的顶点场景往往极端的离谱工程上几乎不可能发生所以解出来往往过于保守。椭球式集合。数学上好看但是会让约束变成二阶锥求解复杂度和实现难度都上了一个档次。多面体集合也叫预算约束集合。在盒式集合的基础上加了一个预算限制即所有不确定参数偏离预测值的总幅度不超过某个值Γ。这个设计很妙因为它允许每个变量偏离但不允许所有变量同时偏离到最极端值Γ越大越保守越小越接近于确定性模型。我们复现的模型使用多面体不确定集合并且每个不确定性变量设置了各自的偏差比例。具体来说风电出力ut的范围是[预测值-0.2×预测值, 预测值0.2×预测值]光伏相同负荷的偏差范围设为±0.1。同时在时域上加了一个总预算约束也就是整个调度周期内每个时段的不确定变量偏离基准值的累计曼哈顿距离被限制住。这样的好处是不会出现“全天所有时段的风光同时剧烈波动”这种不合理场景。为什么这么设计你想想一个系统如果按照“所有不确定性同时拉满”去配置容量那备用成本会高到不可接受。而预算约束集合相当于让不确定性的“总破坏力”保持在一个现实范围内反映的是大数定律和气象条件的平滑效应——实际风光出力虽然单点有波动但不会每个时刻都同时处在极端。3. 从公式到复现模型变换与CCG算法落地全流程3.1 工具选型为什么用PythonGurobi个人观点微电网这种中小规模优化问题PythonGurobi是复现效率最高的组合没有之一。Gurobi的线性规划和对偶单纯形法在高稀疏度矩阵上表现极其稳定而且Python接口写起来和数学表达式几乎一一对应不容易出现“模型翻译错误”。当然有人会问为什么不直接用MatlabYALMIPYALMIP写两阶段鲁棒优化需要手动管理不确定变量和矩阵重构代码长了之后debug的体验很崩溃。Python这边直接用gurobipy加上numpy做矩阵预处理整个代码结构非常清爽。如果是Gurobi不可用的场景用通用的开源求解器替代也是可行的但要特别注意求解时间CCG的主问题是一个混合整数线性规划开源MIP求解器在中等规模下容易卡壳。我们复现的算例规模是24个调度时段2台燃气轮机1台储能1个风电场1个光伏电站1个可调负荷节点。MIP变量约300个约束约500行。这个规模在Gurobi下求解一个主问题大约在0.5到2秒之间整个CCG迭代20次上下总耗时不到1分钟完全在可接受的范围内。3.2 CCG算法的核心流程每一步都在干什么CCG算法的流程可以用五步说清楚这里我结合微电网的应用场景逐步拆解第一步初始化。设置一个初始的不确定场景一般取预测值场景即假设后续时段的风光出力和负荷都不偏离预测值。这样初始化的目的是给算法一个合理的起点。第二步求解主问题MP。把当前已知的所有最坏场景对应的“第二阶段变量”和“第二阶段约束”加入到主问题中。由于每个加入的场景都会生成一组独立的第二阶段变量随着迭代进行主问题规模会逐渐变大。主问题的解得到当前最优的第一阶段决策y*其目标函数值作为下界LB。第三步求解子问题SP。在此不详细展开对偶推导的每个中间步骤但它的核心思路是把第二阶段问题里面的第一层“min”通过KKT条件或对偶理论等价去掉从而得到一个关于不确定性变量u的单层max问题。这个对偶变换是整个算法的数学核心也是最容易出错的地方。常见错误包括内层min问题不是凸的、约束里混入了整数变量导致对偶不正确、对偶乘子的维度对不上约束维度。这里有个实战中的替代思路如果对偶推导属实憋不出来可以用穷举法在不确定集合的顶点集上枚举u的取值然后逐个求解内层的min问题取最大的那个。对于不确定性变量数量不超过10个的算例顶点枚举法虽然笨但结果完全正确适合用来验证对偶推导是否正确。第四步收敛性判断。如果UB和LB之间的相对间隙小于设定的容差比如0.1%算法终止输出当前y*作为最优解。如果不收敛继续第五步。第五步生成新列。取子问题的最优解u*回到第二步把u*对应的第二阶段场景变量和约束加入到主问题中重新求解。这个过程每迭代一次主问题就会多一组变量和约束这也是“列与约束生成”这个名字的来源。为了便于复现我给出主问题部分的关键代码框架这段代码的结构是按Gurobi的Python接口写的# MP主问题第一阶段变量已添加的不确定场景集合 y {} # 第一阶段变量 for t in range(T): y[on_%d % t] mp.addVar(vtypeGRB.BINARY, nameon_%d % t) y[buy_%d % t] mp.addVar(lb-GRB.INFINITY, ubGRB.INFINITY, namebuy_%d % t) x {} # 第二阶段变量按场景k区分索引 for k in range(K): # K是已加入的极端场景数量 for t in range(T): x[pg_%d_%d % (k, t)] mp.addVar(lb0, ubPG_MAX, namepg_%d_%d % (k, t)) # 其他变量类似... # 第二阶段约束需要为每个场景各写一份 for k in range(K): for t in range(T): # 功率平衡约束u_wind[k][t]是场景k在t时段的风电出力数值 mp.addConstr( x[pg_%d_%d % (k, t)] x[pbat_%d_%d % (k, t)] x[buy_%d_%d % (k, t)] u_wind[k][t] u_pv[k][t], GRB.EQUAL, load[t] x[pcur_%d_%d % (k, t)] x[pload_%d_%d % (k, t)] )在实际代码中第二阶段变量通常需要按场景展开场景越多主问题越大这个现象是正常的。3.3 子问题的对偶变换全程手把手解析子问题是整个算法的体力活所在。第二阶段问题的标准形式是min d^T x s.t. Bx ≥ h - Ey - Gu x ≥ 0部分变量可带上下界这个问题的对偶问题写成max μ ≥ 0, μ^T (h - Ey - Gu) s.t. B^T μ ≤ d因为目标函数是μ^T乘以常数向量h - Ey减去μ^T乘以Gu而u本身又是变量所以整个对偶子问题变成一个以μ和u为变量的双线性规划。双线性项来自μ^TGu。这里对偶乘子μ和不确定性变量u相乘使得问题不再是一个标准的线性规划而是一个双线性规划。处理这个双线性项最常用的路子是强对偶理论大M法线性化。因为μ是有界的B^Tμ ≤ du也是有界的u ∈ U是多面体所以可以引入辅助变量z μ·u再用大M变量把双线性约束线性化。具体来说当μ是连续变量、u是连续变量时需要用McCormick包络做松弛当u在顶点取0或1的指示变量时用大M法加Big-M约束即可。我们复现的场景中不确定性变量是连续的风光出力和负荷所以用McCormick包络来处理这部分非线性。这里有一个经验之谈如果对偶子问题直接建模遇到数值稳定性问题Gurobi经常报unsolved或infeasible优先检查的是大M系数是否取得足够大以及不确定集合的边界是否闭合。另外如果一个变量有物理上下界建议直接在变量定义时施加而不是通过约束去限制能减少很多病态。为了减少数值困难我们的做法是给每个不确定性变量的离散取值步长设置一个上限比如风电出力最多取20个离散值这样u变成一个离散集双线性项就转化为离散组合枚举线性化起来要可靠得多。代价是最优值会有微小的离散化误差但对于工程调度而言误差在0.1%以内完全可以接受。3.4 完整迭代框架LB、UB怎么更新主问题解作为LB子问题解作为UB这是CCG最常见也最让人迷糊的地方。我解释一句LB是给所有不确定性场景兜底的成本下界因为主问题只考虑了当前两个场景还没考虑所有可能的坏场景UB是把第一阶段的决策固定下来后在最坏场景下算出的真实成本上限。两者的差越来越小说明模型的“兜底计划”越来越完整。迭代过程中LB和UB的更新方式如下LB -np.inf UB np.inf for it in range(max_iter): # step1: 求解MP mp.optimize() if mp.status ! GRB.OPTIMAL: break LB max(LB, mp.objVal) # step2: 求解SP得到最坏场景u_star和子问题目标obj_sp u_star, obj_sp solve_sp(y_star, uncertain_params) UB min(UB, obj_sp) if (UB - LB) / abs(LB) tol: break # step3: 把u_star作为新场景加入MP add_scenario_to_mp(u_star)注意这里UB的更新用的是min因为每轮求出的最坏场景成本可能比之前低LB用max因为主问题考虑的场景越全下界越接近真实最优值。如果一个代码实现里把UB写成max、把LB写成min那算法永远不收敛这是新手最容易踩的坑之一。4. 复现过程中的严格验证怎么确定你的结果是对的4.1 小规模算例用手推结果锁定正确性刚写完代码千万不要直接上24时段的大算例。我的习惯是先构造一个2时段、单机组、单不确定变量的最小算例用一个极端场景手算预期的调度计划再让程序跑出来比对。举个具体的验证例子假设只有一台燃气轮机2个时段燃气轮机爬坡上限是10MW风电2时段出力分别为5MW和15MW负荷均为10MW机组最大出力20MW接线简单。如果不考虑鲁棒确定性模型会在风电预测为5MW时让机组发5MW在15MW时让机组发0MW。但现在假设风电不确定区间是[预测值-2, 预测值2]则时段1风电最坏情况是3MW需要机组发7MW时段2最坏情况是13MW算下来需要机组发0但由于机组在时段1已经在发7MW因此时段2只需降出力至0这没问题。如果不是2个时段而是连续多个时段爬坡约束就可能把这种简单结论打破你需要专门设计一个能让爬坡约束生效的用例来测试。当这段程序的运行结果和手推完全一致时基本可以确定求解内核没问题再扩展到24时段才有意义。还有一个验证方法是对比确定性模型把不确定集合的波动范围设为0两阶段鲁棒模型退化成普通的确定性经济调度模型此时结果必须和标准的单阶段优化结果完全一致。如果这个一致性都保证不了说明代码里存在约束冲突或变量索引错误。4.2 结果的三层合理性分析算例跑通之后不要急着写总结先对结果做三层合理性分析。第一层看调度计划是否满足物理规律。储能SOC曲线不能出现跳变充放电切换不能过于频繁除非电池退化成本设得极低燃气轮机出力变化率必须在爬坡范围内购售电曲线要和分时电价有明显的“低谷买、高峰卖”关系。如果这些基本规律都违背了先检查约束是否漏加了。第二层看经济性指标是否符合直觉。鲁棒优化的结果成本必然要高于同参数下的确定性模型成本高出多少就是“鲁棒代价”。一般来说当用户批评你的结果太贵时通常是因为不确定集合设置得过大。一个合理的区间是当风电波动范围设为±20%时鲁棒优化成本比确定性模型高5%15%。如果高出30%以上大概率不确定集合的预算参数设过头了。第三层看极限情况下的行为。设计一个极端测试把不确定偏差上限设为0且预算Γ设为0结果应与确定性模型一致。把偏差上限设为80%系统应该出现大量切负荷或购买高昂实时电力的行为。这两个极端行为能帮你确认模型的行为逻辑是健全的。4.3 敏感性分析告诉你为什么要关注Γ而不是盲目调大预算约束Γ是控制不确定集合大小的核心参数复现时建议多做几组Γ的敏感性分析。我从实际实验里拿到的典型结果如下Γ不确定预算鲁棒优化总成本元确定性模型成本元成本增加比例最坏场景切负荷量0即确定性12680126800%0213245126804.5%0513520126806.6%010142861268012.7%1.5%24全天极端158201268024.7%8.2%可以看到Γ从2增加到10成本增加明显但系统还能保证不切负荷到了Γ24相当于全天每个时段都允许极端偏差且同时发生成本暴涨25%还出现了显著切负荷。这说明在实际决策时把Γ设在510之间用6%12%的成本冗余换取全场景可行是工程上可以接受的操作区间。这里多说一句很多论文只报告“最坏情况下的成本”却不想一个更重要的问题实际场景下不确定性没有到达最坏值鲁棒优化的调度方案表现如何如果你感兴趣可以在复现后把CCG求出的第一阶段计划固定住然后用10万个蒙特卡洛抽样场景去验算实际平均成本和最大成本你会发现鲁棒方案的平均成本只比确定性方案高一点点但最坏情况下的成本大幅降低。这就是鲁棒调度的真正价值所在。5. 复现中的常见问题与排查技巧这些坑我替你踩过了5.1 问题速查表按现象找原因我把自己复现过程中遇到的高频问题整理成一张速查表方便你对照排查现象可能原因排查方法主问题无解功率平衡约束写错方向或机组出力上限小于最低负荷去掉鲁棒约束先跑确定性模型验证模型本身是否可行子问题无解或对偶错误第二阶段可行域为空或对偶推导漏掉某个约束检查是否真的对每个约束都写了对应的对偶变量尝试用枚举顶点法验证LB和UB不收敛间隙震荡大M系数设置不当或场景添加逻辑出错检查UB的更新是否存在用了错误的场景索引调大大M系数并加松弛变量结果过于保守成本离谱不确定集合的预算Γ设得过大或偏差比例设得过高减小Γ观察成本变化是否平滑对比确定性模型储能SOC曲线阶梯式跳变储能容量约束没有按时段连续关联或SOC更新公式索引错误打印每个时段的SOC值逐时段检查递推公式求解速度突然变得极慢场景数过多主问题规模膨胀在MP中对已固定的场景参数做预计算删除冗余变量或改用对偶子问题快速求解5.2 最容易出错的三个细节你的代码大概率也栽在这第一个坑是第二阶段变量索引错位。当你把多个最坏场景依次加入主问题时每个场景都对应一套独立的第二阶段变量。如果索引写的不是(k, t)而是(t)后加入的场景会把前面场景的变量覆盖掉结果就乱了。排查方法很简单迭代5次后把主问题里的场景数打印出来看是否严格等于迭代次数加1。第二个坑是储能SOC递推公式在鲁棒场景下的写法。很多复现者习惯把SOC写成单一变量但在CCG里储能每个场景都应有独立的SOC轨迹因为不同场景下充放电功率不同SOC必然不同。主问题的SOC约束需要带上场景索引k。第三个坑是对偶乘子的非负性处理。第二阶段约束既有等式又有不等式等式约束对应的对偶变量是自由变量不等式约束的对偶变量才要求非负。好多复现者图省事把所有对偶变量都加了非负限制结果子问题永远找不到最优解。正确做法是严格区分等式和不等式约束的对偶乘子属性。5.3 微电网经济调度落地的两条实用经验在复现之外我还想分享两条从实际工程项目里得来的经验尤其是矿山微电网这类场景这两条经验能帮你少走很多弯路。第一条矿山微电网的负荷往往是冲击性负荷大型提升机、破碎机启动瞬间的功率冲击可以达到稳态运行的好几倍这在设计不确定集合时必须考虑到。常规微电网负荷不确定集合如果只设±10%的波动范围在矿山场景下是不够的我建议至少按±20%设计并且要把“冲击性负荷发生时段的功率跳变”单独建一个不确定max场景否则鲁棒优化方案在实际运行中依然会频繁越限。第二条风光的预测精度直接决定了不确定集合的设计依据复现时如果有条件可以用历史预测数据和实际数据的误差分布来标定偏差边界而不是拍脑袋定20%。我在做矿山微电网项目时曾用过去两个月的风电和光伏预测误差数据取95%分位数作为偏差上限这样得到的鲁棒方案既不过分保守又留有足够的安全裕度。这种方法比纯靠猜要可靠得多。6. 这类调度系统后续还能怎么扩展任何方法在理论验证完成之后都要接受工程化的考验两阶段鲁棒优化也不例外。基于复现过程中的积累我列几条在自己项目中验证过、性价比比较高的扩展方向第一在多时间尺度上做嵌套。日前用两阶段鲁棒日内滚动用MPC模型预测控制把日前定下来的机组组合作为MPC的边界条件MPC实时修正出力偏差。这种方案的优点是既有鲁棒优化的全局视角又有MPC的实时修正能力实际运行效果非常好。第二把两阶段扩展为多阶段或者引入数据驱动的鲁棒优化。传统鲁棒优化不利用历史数据直接给一个固定集合如果集合偏大就保守偏小就不够鲁棒。数据驱动鲁棒优化通过历史数据构造不确定性的置信集合能有效缓解这个问题。这一方向在研究界和工业界的关注度都在快速增长。第三加入碳交易机制和绿电消纳考核。微电网调度逐渐从“纯成本最小化”走向“成本碳排放绿电消纳等多目标”把碳价、绿证价格作为参数引入目标函数模型形式不用大改但决策逻辑会发生根本变化。以我个人的经验越是复杂的算法落地时越要关注输入数据的质量和不确定集合设计的合理性算法本身反而比较成熟。微电网项目的差异主要在场景适配和数据基础这也是为什么同样的两阶段鲁棒代码在不同项目里效果差异巨大。这套复现方案已经帮你把“算法骨架”搭好了下一步的优化空间在数据和场景这两块而不在代码本身。 两阶段鲁棒优化微电网经济调度：CCG算法复现全解析
1. 项目缘起这不是“套模型”是一场和不确定性较劲的实战做微电网优化调度的人应该都有同感你辛辛苦苦搭了一个确定性经济调度模型把风光预测值往里面一塞算出来一套调度计划以为万事大吉。结果第二天光伏被云遮了半小时风机出力跟预测差了快两成实时平衡市场电价一波动原本算好的“最优”方案瞬间变成了一笔亏钱方案。这种挫败感做过工程的人都能懂。我这次复现的两阶段鲁棒优化经济调度方法解决的就是这个问题。它的核心思路是不在预测值上赌运气而是假设不确定性在某个给定范围内任意变化找一套在所有可能场景下都能保证可行、并且最坏情况下总成本最优的调度方案。说白了就是用一点点经济性代价换取系统运行的安全性边界。这套方法在当前微电网、主动配电网、综合能源系统的研究中是妥妥的热点方向也是工业界从“计划调度”走向“韧性调度”绕不开的一步。这篇博文不是给论文做摘要而是把从数学模型、算法设计、代码复现到结果验证的完整链路记录下来。适合三类人看一是正在复现相关论文、卡在两阶段鲁棒建模细节里的研究生二是想在企业微电网尤其是矿山、工业园区这类对供电可靠性要求苛刻的场景落地鲁棒调度的工程师三是刚开始接触鲁棒优化、想从零理解“列与约束生成”到底在干什么的入门者。如果你是第三类建议把第一节和第二节反复看两遍后面再跟着代码实操。2. 两阶段鲁棒优化到底在干什么先把这个模型掰开了理解2.1 两阶段的结构先定“现在”再应对“未来”两阶段鲁棒优化的名字听起来高大上但拆开看就是一个“现在做决定、未来做调整”的决策框架。拿微电网来举例调度员在日前比如上午10点就需要确定次日的机组启停、与主网的购售电计划这是第一阶段决策特点是决策时不知道确切的分布式电源出力和负荷。等第二天实际运行到某个时刻风光出力、负荷的真实数据暴露了一部分调度员需要在第一阶段决策的基础上对机组出力爬坡、电池充放电功率、切负荷量做实时调整来保证系统功率平衡这是第二阶段决策。这里有几个容易绕晕的点我展开说一下第一阶段决策也叫“here-and-now”决策特点是决策变量在不确定性实现前就必须敲定具有不可更改性。比如火电机组的启停状态你不能说“等风来了再决定启不启机”。第二阶段决策叫“wait-and-see”决策也就是看到不确定性实现值之后再做的适应性调整。它存在的前提是系统有灵活性资源比如储能、可调负荷、机组爬坡能力。鲁棒优化和随机规划最本质的区别就在这里随机规划给不确定性变量指定概率分布然后用期望值目标鲁棒优化不给分布只给定一个不确定集合然后找最坏情况下的最优解。如果你从来没接触过这类模型可以拿买房做类比。第一阶段决策是“签贷款合同”这时你锁定了贷款方案第二阶段决策是“交月供”如果利率上浮不确定性实现你需要调整每月还款结构来保证不违约。鲁棒优化的思路就是提前假设利率可能在某个区间内波动然后找一个无论利率怎么变都不会违约、并且最坏情况下总还款额最低的贷款方案。2.2 模型表达式用数学语言把这个调度问题锁死我们今天复现的微电网拓扑结构是典型的交流微电网带储能系统包含燃气轮机、储能电池、风机、光伏、固定负荷和可调负荷同时和配电网有联络线可以双向购售电。这里我把核心数学模型写出来注意这里用的不是论文里那种缩写极多的套路而是可以直接翻译成代码的写法。目标函数是所有阶段的总成本最小包括第一阶段成本燃气轮机启停成本、日前与主网的购电成本。第二阶段成本燃气轮机燃料成本、某些可调负荷的响应补偿成本、弃风弃光惩罚成本、最坏场景下的切负荷惩罚成本。写成标准的两阶段鲁棒优化形式就是min c1^T y max u∈U min x(z,u) c2^T xs.t. Ay ≥ b 这是第一阶段约束比如机组启停逻辑、最小启停时间约束Bx ≥ d - Cu - Ey 这是第二阶段约束比如功率平衡、机组出力上下限、储能SOC递推公式u ∈ U U就是我们设定的不确定集合。这里必须强调一个理解上的关键卡点整个模型是一个“min-max-min”三层结构。内层的min是第二阶段优化问题max是找到一个最坏的不确定性实现u让总成本最大外层的min是选择第一阶段决策变量y让最坏情况下的成本最小。这个结构直接求解是不可能的因为内层的max会把可行域搞得不可凸所以必须用算法去迭代逼近。目前主流的解法是两个列与约束生成算法也是我们这次复现用的方法。它的思想是先固定一个不确定场景求解主问题得到一个下界然后把主问题的解代入子问题求解一个max-min问题找到真正最坏的不确定性场景再把这个场景作为新的约束列加入到主问题里反复迭代直到上下界收敛。Benders对偶分解思想是把max-min子问题通过对偶变换转化为单层max问题然后用割平面去逼近。这个方法的缺点是当第二阶段是整数变量时不好用。选择CCG而不是Benders原因很直接CCG在处理含整数变量比如可调负荷是否参与响应的0-1状态时不需要做主问题的多次可行割回代收敛速度要快得多。在实际的微电网调度中储能充放电状态、可调负荷的投入状态这类整数变量几乎必然出现所以CCG是更稳妥的选择。2.3 不确定集合的设计决定鲁棒性和经济性的天平很多人在复现论文时把精力全花在求解算法上却忽视了在最前面就要做的一个决定——不确定集合长什么样。不确定集合就是给风光出力和负荷预设一个波动区间形式直接决定了模型的保守程度。常见的有三种盒式集合。最简单每个不确定性变量独立在一个区间内变化。缺点是它的顶点场景往往极端的离谱工程上几乎不可能发生所以解出来往往过于保守。椭球式集合。数学上好看但是会让约束变成二阶锥求解复杂度和实现难度都上了一个档次。多面体集合也叫预算约束集合。在盒式集合的基础上加了一个预算限制即所有不确定参数偏离预测值的总幅度不超过某个值Γ。这个设计很妙因为它允许每个变量偏离但不允许所有变量同时偏离到最极端值Γ越大越保守越小越接近于确定性模型。我们复现的模型使用多面体不确定集合并且每个不确定性变量设置了各自的偏差比例。具体来说风电出力ut的范围是[预测值-0.2×预测值, 预测值0.2×预测值]光伏相同负荷的偏差范围设为±0.1。同时在时域上加了一个总预算约束也就是整个调度周期内每个时段的不确定变量偏离基准值的累计曼哈顿距离被限制住。这样的好处是不会出现“全天所有时段的风光同时剧烈波动”这种不合理场景。为什么这么设计你想想一个系统如果按照“所有不确定性同时拉满”去配置容量那备用成本会高到不可接受。而预算约束集合相当于让不确定性的“总破坏力”保持在一个现实范围内反映的是大数定律和气象条件的平滑效应——实际风光出力虽然单点有波动但不会每个时刻都同时处在极端。3. 从公式到复现模型变换与CCG算法落地全流程3.1 工具选型为什么用PythonGurobi个人观点微电网这种中小规模优化问题PythonGurobi是复现效率最高的组合没有之一。Gurobi的线性规划和对偶单纯形法在高稀疏度矩阵上表现极其稳定而且Python接口写起来和数学表达式几乎一一对应不容易出现“模型翻译错误”。当然有人会问为什么不直接用MatlabYALMIPYALMIP写两阶段鲁棒优化需要手动管理不确定变量和矩阵重构代码长了之后debug的体验很崩溃。Python这边直接用gurobipy加上numpy做矩阵预处理整个代码结构非常清爽。如果是Gurobi不可用的场景用通用的开源求解器替代也是可行的但要特别注意求解时间CCG的主问题是一个混合整数线性规划开源MIP求解器在中等规模下容易卡壳。我们复现的算例规模是24个调度时段2台燃气轮机1台储能1个风电场1个光伏电站1个可调负荷节点。MIP变量约300个约束约500行。这个规模在Gurobi下求解一个主问题大约在0.5到2秒之间整个CCG迭代20次上下总耗时不到1分钟完全在可接受的范围内。3.2 CCG算法的核心流程每一步都在干什么CCG算法的流程可以用五步说清楚这里我结合微电网的应用场景逐步拆解第一步初始化。设置一个初始的不确定场景一般取预测值场景即假设后续时段的风光出力和负荷都不偏离预测值。这样初始化的目的是给算法一个合理的起点。第二步求解主问题MP。把当前已知的所有最坏场景对应的“第二阶段变量”和“第二阶段约束”加入到主问题中。由于每个加入的场景都会生成一组独立的第二阶段变量随着迭代进行主问题规模会逐渐变大。主问题的解得到当前最优的第一阶段决策y*其目标函数值作为下界LB。第三步求解子问题SP。在此不详细展开对偶推导的每个中间步骤但它的核心思路是把第二阶段问题里面的第一层“min”通过KKT条件或对偶理论等价去掉从而得到一个关于不确定性变量u的单层max问题。这个对偶变换是整个算法的数学核心也是最容易出错的地方。常见错误包括内层min问题不是凸的、约束里混入了整数变量导致对偶不正确、对偶乘子的维度对不上约束维度。这里有个实战中的替代思路如果对偶推导属实憋不出来可以用穷举法在不确定集合的顶点集上枚举u的取值然后逐个求解内层的min问题取最大的那个。对于不确定性变量数量不超过10个的算例顶点枚举法虽然笨但结果完全正确适合用来验证对偶推导是否正确。第四步收敛性判断。如果UB和LB之间的相对间隙小于设定的容差比如0.1%算法终止输出当前y*作为最优解。如果不收敛继续第五步。第五步生成新列。取子问题的最优解u*回到第二步把u*对应的第二阶段场景变量和约束加入到主问题中重新求解。这个过程每迭代一次主问题就会多一组变量和约束这也是“列与约束生成”这个名字的来源。为了便于复现我给出主问题部分的关键代码框架这段代码的结构是按Gurobi的Python接口写的# MP主问题第一阶段变量已添加的不确定场景集合 y {} # 第一阶段变量 for t in range(T): y[on_%d % t] mp.addVar(vtypeGRB.BINARY, nameon_%d % t) y[buy_%d % t] mp.addVar(lb-GRB.INFINITY, ubGRB.INFINITY, namebuy_%d % t) x {} # 第二阶段变量按场景k区分索引 for k in range(K): # K是已加入的极端场景数量 for t in range(T): x[pg_%d_%d % (k, t)] mp.addVar(lb0, ubPG_MAX, namepg_%d_%d % (k, t)) # 其他变量类似... # 第二阶段约束需要为每个场景各写一份 for k in range(K): for t in range(T): # 功率平衡约束u_wind[k][t]是场景k在t时段的风电出力数值 mp.addConstr( x[pg_%d_%d % (k, t)] x[pbat_%d_%d % (k, t)] x[buy_%d_%d % (k, t)] u_wind[k][t] u_pv[k][t], GRB.EQUAL, load[t] x[pcur_%d_%d % (k, t)] x[pload_%d_%d % (k, t)] )在实际代码中第二阶段变量通常需要按场景展开场景越多主问题越大这个现象是正常的。3.3 子问题的对偶变换全程手把手解析子问题是整个算法的体力活所在。第二阶段问题的标准形式是min d^T x s.t. Bx ≥ h - Ey - Gu x ≥ 0部分变量可带上下界这个问题的对偶问题写成max μ ≥ 0, μ^T (h - Ey - Gu) s.t. B^T μ ≤ d因为目标函数是μ^T乘以常数向量h - Ey减去μ^T乘以Gu而u本身又是变量所以整个对偶子问题变成一个以μ和u为变量的双线性规划。双线性项来自μ^TGu。这里对偶乘子μ和不确定性变量u相乘使得问题不再是一个标准的线性规划而是一个双线性规划。处理这个双线性项最常用的路子是强对偶理论大M法线性化。因为μ是有界的B^Tμ ≤ du也是有界的u ∈ U是多面体所以可以引入辅助变量z μ·u再用大M变量把双线性约束线性化。具体来说当μ是连续变量、u是连续变量时需要用McCormick包络做松弛当u在顶点取0或1的指示变量时用大M法加Big-M约束即可。我们复现的场景中不确定性变量是连续的风光出力和负荷所以用McCormick包络来处理这部分非线性。这里有一个经验之谈如果对偶子问题直接建模遇到数值稳定性问题Gurobi经常报unsolved或infeasible优先检查的是大M系数是否取得足够大以及不确定集合的边界是否闭合。另外如果一个变量有物理上下界建议直接在变量定义时施加而不是通过约束去限制能减少很多病态。为了减少数值困难我们的做法是给每个不确定性变量的离散取值步长设置一个上限比如风电出力最多取20个离散值这样u变成一个离散集双线性项就转化为离散组合枚举线性化起来要可靠得多。代价是最优值会有微小的离散化误差但对于工程调度而言误差在0.1%以内完全可以接受。3.4 完整迭代框架LB、UB怎么更新主问题解作为LB子问题解作为UB这是CCG最常见也最让人迷糊的地方。我解释一句LB是给所有不确定性场景兜底的成本下界因为主问题只考虑了当前两个场景还没考虑所有可能的坏场景UB是把第一阶段的决策固定下来后在最坏场景下算出的真实成本上限。两者的差越来越小说明模型的“兜底计划”越来越完整。迭代过程中LB和UB的更新方式如下LB -np.inf UB np.inf for it in range(max_iter): # step1: 求解MP mp.optimize() if mp.status ! GRB.OPTIMAL: break LB max(LB, mp.objVal) # step2: 求解SP得到最坏场景u_star和子问题目标obj_sp u_star, obj_sp solve_sp(y_star, uncertain_params) UB min(UB, obj_sp) if (UB - LB) / abs(LB) tol: break # step3: 把u_star作为新场景加入MP add_scenario_to_mp(u_star)注意这里UB的更新用的是min因为每轮求出的最坏场景成本可能比之前低LB用max因为主问题考虑的场景越全下界越接近真实最优值。如果一个代码实现里把UB写成max、把LB写成min那算法永远不收敛这是新手最容易踩的坑之一。4. 复现过程中的严格验证怎么确定你的结果是对的4.1 小规模算例用手推结果锁定正确性刚写完代码千万不要直接上24时段的大算例。我的习惯是先构造一个2时段、单机组、单不确定变量的最小算例用一个极端场景手算预期的调度计划再让程序跑出来比对。举个具体的验证例子假设只有一台燃气轮机2个时段燃气轮机爬坡上限是10MW风电2时段出力分别为5MW和15MW负荷均为10MW机组最大出力20MW接线简单。如果不考虑鲁棒确定性模型会在风电预测为5MW时让机组发5MW在15MW时让机组发0MW。但现在假设风电不确定区间是[预测值-2, 预测值2]则时段1风电最坏情况是3MW需要机组发7MW时段2最坏情况是13MW算下来需要机组发0但由于机组在时段1已经在发7MW因此时段2只需降出力至0这没问题。如果不是2个时段而是连续多个时段爬坡约束就可能把这种简单结论打破你需要专门设计一个能让爬坡约束生效的用例来测试。当这段程序的运行结果和手推完全一致时基本可以确定求解内核没问题再扩展到24时段才有意义。还有一个验证方法是对比确定性模型把不确定集合的波动范围设为0两阶段鲁棒模型退化成普通的确定性经济调度模型此时结果必须和标准的单阶段优化结果完全一致。如果这个一致性都保证不了说明代码里存在约束冲突或变量索引错误。4.2 结果的三层合理性分析算例跑通之后不要急着写总结先对结果做三层合理性分析。第一层看调度计划是否满足物理规律。储能SOC曲线不能出现跳变充放电切换不能过于频繁除非电池退化成本设得极低燃气轮机出力变化率必须在爬坡范围内购售电曲线要和分时电价有明显的“低谷买、高峰卖”关系。如果这些基本规律都违背了先检查约束是否漏加了。第二层看经济性指标是否符合直觉。鲁棒优化的结果成本必然要高于同参数下的确定性模型成本高出多少就是“鲁棒代价”。一般来说当用户批评你的结果太贵时通常是因为不确定集合设置得过大。一个合理的区间是当风电波动范围设为±20%时鲁棒优化成本比确定性模型高5%15%。如果高出30%以上大概率不确定集合的预算参数设过头了。第三层看极限情况下的行为。设计一个极端测试把不确定偏差上限设为0且预算Γ设为0结果应与确定性模型一致。把偏差上限设为80%系统应该出现大量切负荷或购买高昂实时电力的行为。这两个极端行为能帮你确认模型的行为逻辑是健全的。4.3 敏感性分析告诉你为什么要关注Γ而不是盲目调大预算约束Γ是控制不确定集合大小的核心参数复现时建议多做几组Γ的敏感性分析。我从实际实验里拿到的典型结果如下Γ不确定预算鲁棒优化总成本元确定性模型成本元成本增加比例最坏场景切负荷量0即确定性12680126800%0213245126804.5%0513520126806.6%010142861268012.7%1.5%24全天极端158201268024.7%8.2%可以看到Γ从2增加到10成本增加明显但系统还能保证不切负荷到了Γ24相当于全天每个时段都允许极端偏差且同时发生成本暴涨25%还出现了显著切负荷。这说明在实际决策时把Γ设在510之间用6%12%的成本冗余换取全场景可行是工程上可以接受的操作区间。这里多说一句很多论文只报告“最坏情况下的成本”却不想一个更重要的问题实际场景下不确定性没有到达最坏值鲁棒优化的调度方案表现如何如果你感兴趣可以在复现后把CCG求出的第一阶段计划固定住然后用10万个蒙特卡洛抽样场景去验算实际平均成本和最大成本你会发现鲁棒方案的平均成本只比确定性方案高一点点但最坏情况下的成本大幅降低。这就是鲁棒调度的真正价值所在。5. 复现中的常见问题与排查技巧这些坑我替你踩过了5.1 问题速查表按现象找原因我把自己复现过程中遇到的高频问题整理成一张速查表方便你对照排查现象可能原因排查方法主问题无解功率平衡约束写错方向或机组出力上限小于最低负荷去掉鲁棒约束先跑确定性模型验证模型本身是否可行子问题无解或对偶错误第二阶段可行域为空或对偶推导漏掉某个约束检查是否真的对每个约束都写了对应的对偶变量尝试用枚举顶点法验证LB和UB不收敛间隙震荡大M系数设置不当或场景添加逻辑出错检查UB的更新是否存在用了错误的场景索引调大大M系数并加松弛变量结果过于保守成本离谱不确定集合的预算Γ设得过大或偏差比例设得过高减小Γ观察成本变化是否平滑对比确定性模型储能SOC曲线阶梯式跳变储能容量约束没有按时段连续关联或SOC更新公式索引错误打印每个时段的SOC值逐时段检查递推公式求解速度突然变得极慢场景数过多主问题规模膨胀在MP中对已固定的场景参数做预计算删除冗余变量或改用对偶子问题快速求解5.2 最容易出错的三个细节你的代码大概率也栽在这第一个坑是第二阶段变量索引错位。当你把多个最坏场景依次加入主问题时每个场景都对应一套独立的第二阶段变量。如果索引写的不是(k, t)而是(t)后加入的场景会把前面场景的变量覆盖掉结果就乱了。排查方法很简单迭代5次后把主问题里的场景数打印出来看是否严格等于迭代次数加1。第二个坑是储能SOC递推公式在鲁棒场景下的写法。很多复现者习惯把SOC写成单一变量但在CCG里储能每个场景都应有独立的SOC轨迹因为不同场景下充放电功率不同SOC必然不同。主问题的SOC约束需要带上场景索引k。第三个坑是对偶乘子的非负性处理。第二阶段约束既有等式又有不等式等式约束对应的对偶变量是自由变量不等式约束的对偶变量才要求非负。好多复现者图省事把所有对偶变量都加了非负限制结果子问题永远找不到最优解。正确做法是严格区分等式和不等式约束的对偶乘子属性。5.3 微电网经济调度落地的两条实用经验在复现之外我还想分享两条从实际工程项目里得来的经验尤其是矿山微电网这类场景这两条经验能帮你少走很多弯路。第一条矿山微电网的负荷往往是冲击性负荷大型提升机、破碎机启动瞬间的功率冲击可以达到稳态运行的好几倍这在设计不确定集合时必须考虑到。常规微电网负荷不确定集合如果只设±10%的波动范围在矿山场景下是不够的我建议至少按±20%设计并且要把“冲击性负荷发生时段的功率跳变”单独建一个不确定max场景否则鲁棒优化方案在实际运行中依然会频繁越限。第二条风光的预测精度直接决定了不确定集合的设计依据复现时如果有条件可以用历史预测数据和实际数据的误差分布来标定偏差边界而不是拍脑袋定20%。我在做矿山微电网项目时曾用过去两个月的风电和光伏预测误差数据取95%分位数作为偏差上限这样得到的鲁棒方案既不过分保守又留有足够的安全裕度。这种方法比纯靠猜要可靠得多。6. 这类调度系统后续还能怎么扩展任何方法在理论验证完成之后都要接受工程化的考验两阶段鲁棒优化也不例外。基于复现过程中的积累我列几条在自己项目中验证过、性价比比较高的扩展方向第一在多时间尺度上做嵌套。日前用两阶段鲁棒日内滚动用MPC模型预测控制把日前定下来的机组组合作为MPC的边界条件MPC实时修正出力偏差。这种方案的优点是既有鲁棒优化的全局视角又有MPC的实时修正能力实际运行效果非常好。第二把两阶段扩展为多阶段或者引入数据驱动的鲁棒优化。传统鲁棒优化不利用历史数据直接给一个固定集合如果集合偏大就保守偏小就不够鲁棒。数据驱动鲁棒优化通过历史数据构造不确定性的置信集合能有效缓解这个问题。这一方向在研究界和工业界的关注度都在快速增长。第三加入碳交易机制和绿电消纳考核。微电网调度逐渐从“纯成本最小化”走向“成本碳排放绿电消纳等多目标”把碳价、绿证价格作为参数引入目标函数模型形式不用大改但决策逻辑会发生根本变化。以我个人的经验越是复杂的算法落地时越要关注输入数据的质量和不确定集合设计的合理性算法本身反而比较成熟。微电网项目的差异主要在场景适配和数据基础这也是为什么同样的两阶段鲁棒代码在不同项目里效果差异巨大。这套复现方案已经帮你把“算法骨架”搭好了下一步的优化空间在数据和场景这两块而不在代码本身。