数据中心微网两阶段鲁棒规划：CCG算法Matlab复现与实操

刚开始接触“考虑灵活性的数据中心微网两阶段鲁棒规划”这个课题时，我第一反应是头大。这几乎集齐了电力系统优化领域最硬核的几个点：微网容量规划、数据中心柔性负荷建模、鲁棒优化下的min-max-min结构，还得用Matlab把整套CCG（列与约束生成）迭代框架搭出来。不少做这个方向的朋友都卡在半路，要么模型推倒重来，要么算法在某个循环里陷入死锁。所以这篇文章我不打算空谈理论，而是把自己复现EI论文全过程的核心细节、踩坑记录和关键代码逻辑完整梳理一遍，希望能帮同样在死磕这个方向的同学节省一周时间。

这个项目解决的问题非常具体：数据中心微网在做容量规划（装多少光伏、储能、柴油机）时，必须同时考虑两个层面的不确定性——光伏出力的随机波动和数据中心负荷的动态变化。传统确定性规划只取一个场景，容易在极端天气或负荷突增时“翻车”；两阶段鲁棒规划则通过构造最恶劣场景，确保任何不确定区间内的波动都不会导致系统失负荷。换句通俗的话说，第一阶段先把该建的东西定下来（投资决策），第二阶段等不确定“牌”亮出来后，再安排储能充放电、柴油机启停和柔性负荷调动这些运行手段来兜底。这跟我们买房做预算一样，先确定首付和贷款上限（投资约束），再根据每个月的收入波动预留应急资金（运行调整），只不过这里的“最坏情况”由不确定性集合精确定义。

对于刚开始入门的同学，我的建议是先别急着碰代码，把这个项目拆成三块：数据中心柔性资源的物理建模、两阶段鲁棒规划的数学原理、以及CCG算法在Matlab中的落地实现。三者缺一不可。我从复现的角度出发，按照“思路拆解—数学建模—代码实现—问题排查”的顺序逐步展开，中间包含大量可直接套用的代码片段和参数配置方案。

1. 项目整体设计与核心思路

1.1 数据中心微网的特殊性在哪

传统微网规划研究得已经很透了，但数据中心微网的难点在于负荷侧独特。数据中心本身有大量机柜、IT设备和冷却系统，这些负荷并非一成不变。通过服务器负载迁移、虚拟化调度、以及冷却系统温度的适度调整，可以实现5%-20%的柔性调节能力。这样一来，微网规划就不能再把负荷当作刚性值来处理，而是要把“可调范围”纳入不确定性框架。这让问题从单纯的“电源侧配置”升级成了“源-荷”双侧协同规划，同时还要考虑储能和柴油机作为紧急备用。

举例来说，一台满载的服务器功耗可能是额定值的90%，通过负载均衡把部分任务切到其他时段或空闲节点，就能腾出一部分可削减负荷。数据中心管理者甚至可以通过延迟非实时任务（比如批量数据分析、模型训练）的执行时间，换取用电成本的下降。这些灵活性在规划阶段如果完全不考虑，就会导致光伏配得过大、储能容量虚高，整体一次性投资增加10%-30%。实际论文中通常用一个“灵活调节因子”描述负荷可调上限，并结合数据中心业务允许的延迟时长来设定时移约束。

1.2 为什么选择两阶段鲁棒而非随机规划

在做不确定性规划时，最常见的三类框架是随机规划、鲁棒优化和分布鲁棒。当时我看了那篇EI论文的摘要后，很快理解了作者选择两阶段鲁棒的理由。随机规划需要为光伏和负荷都假设出精确的概率密度函数（比如Beta分布、正态分布），然后通过蒙特卡洛抽样生成大量场景，计算代价高，而且一旦假设分布与实际偏差过大，规划结果的可信度就打折扣。鲁棒优化则放弃精确分布，只要求不确定参数的波动区间，对应到工程上就是“光伏出力最低可能到多少”“负荷最高可能涨到多少”，这种数据在规划阶段其实是比较容易获取的。

两阶段鲁棒的核心逻辑可以直观地概括为：先在“今天”做投资决策，再假想“明天”自然选择了最不利的光伏和负荷组合，此时做运行调整仍能保证不失负荷。用数学语言描述就是min-max-min结构，即最小化投资成本与最恶劣场景下运行成本之和，其中外层min对应第一阶段的决策变量，max对应不确定性集合内寻找最恶劣场景，内层min则是第二阶段运行变量的调整。这个结构的好处在于，它天然嵌入了“适应性”特征——不像单阶段鲁棒那样在同一个统一的决策动作中被动地接受全部不确定性，第二阶段可以在观察到不确定性的具体值后就地进行响应。

在Matlab代码实现上，两阶段鲁棒通常需要转化为可求解的主问题-子问题迭代。第一部分的主问题是一个包含“割平面”约束的混合整数线性规划，用来确定第一阶段决策变量和成本下界；第二部分的子问题则用来在固定的第一阶段决策下寻找最恶劣场景，并返回成本上界。若上下界间隙小于设定容差（比如1%），则收敛；否则把子问题中识别出的最恶劣场景作为新的场景变量加进主问题，继续迭代。这就是CCG列与约束生成算法的核心流程。

2. 数学建模与原理拆解

2.1 第一阶段投资决策变量的设计

规划层面的决策变量相对直观，主要包括：光伏装机容量（kW）、储能额定功率（kW）、储能额定容量（kWh）、柴油发电机台数与容量（kW），以及各设备是否建设的0-1状态变量。目标函数的第一项就是这些设备的年化投资成本，需要考虑设备寿命、折现率以及单位造价。比如光伏单位造价大概在3500-4500元/kW，磷酸铁锂储能系统的单位容量成本目前在900-1300元/kWh，单位功率成本约1000-1500元/kW，柴油发电机的成本相对便宜，但会带来燃料和维护费用。

在建模时需要特别注意的一点是，决策变量之间的耦合约束非常多。比如储能功率与容量并非独立，通常要满足最大充放电功率为额定功率倍数（例如0.25C）的约束；光伏的占地面积、楼顶可用面积也是隐含约束。更关键的是，如果储能配置偏小，第二阶段运行时就没有足够的调节空间来应对极端场景，会导致规划结果过于乐观，后评估阶段发现某些场景失负荷严重。因此编程序时第一阶段变量不仅要写进目标函数，还要以参数形式出现在第二阶段合法性约束中，这就要求代码中频繁使用M（大M）法构建逻辑关联条件。

在YALMIP环境中定义第一阶段决策变量，我一般建议区分binvar和sdpvar两种类型，规划变量的维度和索引在仿真前需要与时间序列长度（比如8760小时或典型日场景）严格对齐。这里有个容易出错的点：在YALMIP里，如果你把同一个变量既放在投资目标里，又延续到第二阶段约束中，要确保变量类别的统一，避免出现维度不匹配报错。

2.2 第二阶段运行模型与min-max-min框架

第二阶段运行层是在投资决策给定的前提下，求解各时段内各设备的出力计划，目标函数是运行成本最小化，包括从电网购电成本、柴油机燃料成本、以及储能充放电折旧成本。约束条件涵盖节点功率平衡、储能SOC递推方程、光伏出力上限、输电线路容量、数据中心柔性负荷可调范围等。

在子问题里寻找“最恶劣场景”时，我们需要同时选取光伏出力和数据中心负荷的极端组合。直观上是最低光伏出力叠加最高负荷。但由于网络中可能存在多节点、多时段的复杂交互，这个“组合”不一定只在区间端点取得，所以必须借助对偶理论把内层min替换成对偶形式，把双层问题化简为单层最大化问题。经过对偶变换后的子问题本质上是一个双线性项（不确定性变量乘以对偶乘子）构成的非凸问题，工程上通常采用线性化技巧或者场景枚举的方法处理。CCG算法的高明之处在于，它不需要显式枚举所有组合，而是通过迭代只引入真正起作用的“极点”场景，从而显著缩减计算量。

“这里补充一个我对两阶段鲁棒的理解：两阶段鲁棒的灵活性来自第二阶段可以‘看牌出牌’。如果模型退化成单阶段，那就意味着所有决策在不确定性实现前一次性定死，相当于蒙着眼下棋，结果自然保守得多。两阶段结构引入的正是这种‘自适应’能力，它也是数据中心微网柔性负荷能发挥价值的关键前提。”

2.3 不确定性集合的构建方式

电能系统鲁棒优化中最常用的不确定性集合是盒式区间集合，即PV出力在预测值上下一定百分比范围波动，负荷在预测值上下一定百分比浮动。单纯盒式集合的缺点是所有不确定参数可以同时取到最坏值，所得到的规划结果偏保守。为了在鲁棒性和经济性之间取得折中，工程上一般引入预算约束Γ，这也是论文中最常见的处理方式。Γ表示单位时间内不确定参数与预测值偏差绝对值的总和上限，比如一天24个时段，如果Γ取12，则意味着最多只有12个时段能同时达到极端场景；这个思想由Bertsimas和Sim提出，是鲁棒优化领域的标准做法。

在Matlab中构建不确定性集合时，我一般实现以下两类约束。第一类是对每个时段各不确定参数单独加界限约束；第二类是在时间维度上加入总和预算约束，本质上是一个L1范数约束。这个约束会进入max部分，需要对偶化处理。若不确定参数是连续区间，对偶后会出现一组新的非负乘子变量，代码中需要为它们显式分配变量内存。实际测试下来，Γ取全天时段数的50%-70%时结果误差可接受，且不会让运算量显著增大。

3. Matlab代码实现与实操过程

3.1 工具箱选型与求解环境配置

我使用的是YALMIP R20230922版本配合Cplex 12.9，在Matlab R2022a上运行。坦白说，如果你一开始成本预算紧张，用Gurobi也行，两者都能求解混合整数规划问题。但做鲁棒优化子问题时，YALMIP有较方便的双线性项处理机制，配合yalmiptime和implies等命令可以极大减少建模工作量。此外，求解器的选择至关重要，特别是子问题中大量0-1变量和连续变量混合时，Cplex在分支切割算法上的表现通常比开源求解器好，否则很容易在迭代后期出现计算时间爆炸。

环境配置上，需要确认YALMIP的路径已经正确加入Matlab工作区；同时需要确保Cplex的Java接口版本和Matlab的Java版本兼容，这个问题在新版本Matlab（R2022b以上）安装旧版Cplex时尤其常见。建议先跑通自带的demo，比如yalmiptest，确认MILP求解器能被正常识别，再进入正式建模环节。

3.2 代码框架与CCG迭代流程

整体代码结构我分为四个模块。第一模块是参数设置，包括典型的日数据、设备参数、单位成本等。第二模块是在YALMIP中定义主问题（MP），初始只需要包含第一阶段决策变量和一些基本约束；在这个环节需要额外加上一个dummy场景（比如各不确定参数都取预测值），否则主问题约束数量不足，CCG第一轮迭代无法启动。第三模块是子问题（SP），在给定第一阶段解后，固定x，考虑不确定变量的最大化问题和运行层的最小化问题耦合求解。第四模块是外层循环控制，负责更新上界和下界并判断收敛是否达标。

主问题的展开形式有一个很关键的细节，每轮迭代引入一个新的场景k，同时新增一组辅助变量和第二阶段的运行变量。在新场景中，第二阶段变量需要重新定义，但共享第一阶段变量。这个结构很容易写错，我一开始就是在这里搞混了变量索引，导致每轮迭代主问题的变量数量没有增加，但约束数量却在飙升，最后出现数值怪异。实际上，每轮新增场景后，第二阶段变量的数量应该同步增加，这样才能保证新场景对应独立的运行策略。

子问题的处理核心是把“min-max”转换成“max”形式。常用的技术是先保留内层min问题作为带参线性规划，利用对偶理论获得其对偶结构，然后与max问题合并。对偶变换过程中，需要十分注意等式约束的对偶变量（乘子）没有非负要求，不等式约束的乘子必须非负，否则推导出来的KKT条件和对偶表达式会彻底错掉。为了验证对偶是否正确，我曾用一个小型随机场景将子问题直接枚举求解（离散几个不确定水平）和对偶求解的结果做过对比，确保间隙控制在10的负6次方量级以内再加进CCG框架。

3.3 核心代码逻辑与参数配置建议

下面是主问题的YALMIP核心逻辑框架（约简版），读者可借此对照理解整体结构。

% 主问题初始化阶段，定义了建设变量与运行变量的基本关系 x_inv = binvar(N_device, 1); % 设备建设0-1变量 cap = sdpvar(N_device, 1); % 连续容量变量，MW Constraints = []; obj = sum(Cap_cost .* cap + Fix_cost .* x_inv); % 每轮CCG更新后追加场景k相关变量 for k = 1:K p_stor{k} = sdpvar(T, 1); % 储能出力 soc{k} = sdpvar(T, 1); % 荷电状态 x_curt{k} = sdpvar(T, 1); % 弃光量 x_ls{k} = sdpvar(T, 1); % 负荷削减量 Constraints = [Constraints, ... power_balance_k, soc_k_constraint, ... capacity_link_with_inv(x_inv,cap,k)]; obj = obj + sum(oper_cost_k(:) .* ... [p_stor{k}; x_curt{k}; x_ls{k}]); end optimize(Constraints, obj, sdpsettings('solver','cplex','verbose',2));

子问题的关键思路我在代码中的实现方式如下：把投资收益部分的常数项剥离后，先把运行层问题格式化为f'*y，约束写作A_eq*y = b_eq + F*z_unc的形式。其中z_unc是不确定向量（光伏出力与负荷需求合并）。由于A_eq、F都不包含不确定变量，对偶化后的约束右侧会自然出现A_eq'*lambda = f以及F'*lambda + M_aux等结构表达式。利用YALMIP的dual或kkt命令可以辅助验证，但真心建议不要过分依赖自动工具，手动写下推导过程比出错后排查快得多。

关于不确定集合的预算约束Γ，我用以下约束实现：

Constraints = [Constraints, sum(abs(z_unc - z_nominal) ./ z_range) <= Gamma];

这个约束在子问题的对偶化处理中会产生额外乘子，需要谨慎处理。如果直接用YALMIP进行求解，部分非线性约束可能导致求解失败，我建议先利用显式线性变换把绝对值约束展开为两组不等式，再进入模型。

3.4 典型配置参数与算例设计

下面分享一组我在复现中实测可用的参数配置，供参考。典型日选取一个包含夏季高负荷和冬季低光伏的极端场景，时间分辨率取1小时，一个周期24小时。

负荷基线我用的是实际调研数据，数据中心峰值负荷在2.8MW左右，可调范围设定为15%。光伏在夏季午间出力可达到额定容量的85%，冬季阴天则只有20%。储能的充放电效率取95%，SOC运行范围控制在10%-90%。这些数值看似简单，但对结果影响极大，比如储能效率从95%改为90%，最优配置可能就会有明显变化。

4. 常见问题与排查技巧实录

4.1 CCG迭代不收敛或收敛极慢

这是我复现过程中花时间最多的一部分。CCG理论上在有限步内收敛，但实际可能出现振荡。常见原因是主问题里新场景加入后，没有同时引入足够的辅助变量，导致新场景对应的运行变量复用上一轮旧值，产生数值死锁。解决办法是每轮迭代严格检查新增变量个数，并在循环开始时使用clear指令清理YALMIP内部的缓存变量。

另一个容易忽略的问题是容差设置。我对CCG迭代一般设置tol = 1e-3，也就是上下界相对间隙小于千分之一时终止。如果一直不收敛，可以把容差放宽到1e-2观察结果稳定性；如果仍在振荡，则需要回头检查子问题是否能正常求解出小于零的目标值（当原问题不可行时，对偶问题可能无界或不可行）。

4.2 子问题对偶变换推导错误

这是建模环节最隐蔽的错误源。我在推导子问题时，一开始把光伏出力的约束写反方向，导致对偶乘子符号错误，最终识别出的“最恶劣场景”变成了“最优场景”，CCG上下界出现倒挂。调试这个bug时，我采用了场景枚举法做交叉验证：将不确定性集合离散为5个等级，直接遍历全部组合求解原问题，再与对偶化后的子问题结果对比。两者结果一致后，才放心继续往后推进。

还有一点，处理盒式不确定性集合时，最恶劣场景通常出现在顶点，所以枚举端点且保证每个端点都合法性，是检查对偶模型的有效手段。不要只停留在单个场景，最好对多个固定第一阶段决策反复验证。

4.3 YALMIP建模中的数值警告与求解器异常

在运行过程中，我频繁遇到“NaN or Inf”警报，根因往往是大M值设置不合理。把M值设得过大（比如1e6）会导致混合整数规划求解时出现数值病态；设得过小则会错误削减可行域。我的经验是在所有主从耦合约束中，M值取对应容量上限的2-3倍即可。比如储能出力上限是500kW，M取1000-1500就足够。

另外，由于数据中心微网规划模型内包含大量二进制变量（每个待选站址、每个候选设备类型），MILP求解规模会非常大。强烈建议启用求解器的分支预削减策略，并在求解前对变量顺序做调整，把优先关注的核心容量变量排在前面。这可以显著减少分支搜索树宽度，实测下来求解时间能缩短30%-50%。

5. 复现延展方向

5.1 引入分布鲁棒优化框架

在把两阶段鲁棒跑通之后，可以考虑升级到分布鲁棒优化。它介于随机规划与鲁棒优化之间，只需要利用历史数据构造模糊集，不需要给出精确概率分布。这样做的优势是规划结果对分布误差的敏感度更低。具体做法通常是定义基于Wasserstein距离的模糊集，即所有与经验分布距离在一定半径内的分布都纳入考虑，然后在这个集合中搜索最恶劣分布的期望成本。用CCG处理分布鲁棒优化问题时，子问题的内层决策变量会从场景中选择分布，这比普通两阶段鲁棒多一层对偶变量，但整体框架仍然兼容。

从复现到改进这一步是很多论文创新的基础。比如我这版代码里，如果稍作修改，把盒式不确定性集合替换成基于Wasserstein球的模糊集，再在子问题中额外加一个分布参数变量，就能得到分布鲁棒版本的结果。原来两阶段鲁棒的最恶劣场景是单个点，现在的最恶劣场景可能是一个分布，规划外层配置会更具稳健性。

5.2 接入多场景时序数据建模

另一个值得尝试的方向是把单一典型日扩展到多季节典型日集。EI论文里通常用聚类算法（K-means或谱聚类）从全年365天数据中筛出4-6个代表性场景，每个场景包含24小时的光伏和负荷曲线，并附上权重。在此基础上做带场景权重的两阶段随机鲁棒混合模型，规划精度会比单日模型好很多。实现时需要注意不同场景间第一阶段变量完全共享，但第二阶段变量各自独立，CCG迭代的形式框架仍然不变。

5.3 加入储能寿命损耗的精细化建模

数据中心微网规划对储能依赖度很高，而储能寿命损耗与充放电深度、循环次数强相关。在两阶段鲁棒中，如果把储能寿命损耗建模为深度放电的非线性函数，就会破坏线性规划结构。工程上常用的折中方案是引入“等效循环寿命”约束——每进行一次深度充放电折算为若干标准循环，并把累计循环数作为规划周期内的上限约束。这一约束可以用分段线性函数逼近，从而保留整个问题为MILP的可解性。

6. 实操心得与收尾

整个复现过程下来，我最大的体会是：两阶段鲁棒规划本身不是一个个孤立的算法，而是一套“建模-变换-求解”的完整方法论。很多初学者一上来就在Matlab中堆代码，遇到错误一头雾水。实际上，在写第一行代码之前，把主问题和子问题的数学表达式都手动推导一遍，包括每一处对偶变换，是非常必要的。我在推导“子问题对偶”这一步反复验算，还对输出进行了交叉校验，这些前期工作省下了后面几天的调试时间。

最后再分享一个小技巧：给CCG循环中每一轮迭代的上下界值打点记录并绘制收敛曲线。这样不仅能在投稿论文里作为收敛性佐证，还能在算法振荡时快速定位是哪一轮出了问题。这个习惯我一直沿用到现在，很多优化类项目都能用得上。如果后续还要处理更复杂的不确定性模型，这套“确定性等效+迭代割平面”的基本功同样管用。