风光互补制氢合成氨这个方向，近两年在能源系统优化里热度一直很高。我前阵子刚把一个并网/离网双模式下的容量-调度联合优化模型完整复现了一遍，求解用的是Matlab调用Cplex，跑通之后回头看看，整个模型的搭建思路、约束处理方式、还有求解器调参的细节，都比预想中要讲究。这篇文章就把整个过程拆开讲清楚，从物理系统建模到Cplex落地求解，再到并网离网结果对比，适合正在做新能源制氢、电-氢-氨耦合系统规划课题的研究生，也适合想用Matlab+Cplex做优化调度的工程师参考。
1. 这个优化问题到底在求什么
1.1 系统的物理架构与能量流
先把这个系统的家底摸清楚。风光互补制氢合成氨系统，物理上包含五个核心单元：风电机组、光伏阵列、电解水制氢装置、储氢罐、氨合成塔。在并网模式下，还会和外部电网发生功率交换，所以电网联络线也算一个边界环节。
能量流的方向其实很清晰：风、光发电给电解槽供电，电解槽产出的氢气进入储氢罐缓冲，储氢罐出来的氢气与空分装置提供的氮气按3:1比例进入氨合成塔，生成氨产品。中间每一步都有能量损失和物料损耗，这些细节最终都会转化为模型里的效率系数和物料平衡约束。
系统调制的核心矛盾在于时间尺度的不匹配。风光出力的波动是分钟到小时级别的，电解槽动态响应相对快，但氨合成塔属于化工装置，要求进料流量平稳，不能频繁启停，也不适合深度变负荷。储氢罐就是用来解耦电解制氢和氨合成这两个动态特性差异巨大的环节。
1.2 容量优化和调度优化为什么要联立
这个问题拆成两层来看会更清楚。
容量优化解决的是"一次投多少钱、装多大设备"的问题，决策变量包括风机装机台数、光伏装机容量、电解槽额定功率、储氢罐有效容积、氨合成塔设计产能，以及离网模式下需要配套的储能电池容量。这些变量决定系统的初始投资，属于长期决策。
调度优化解决的是"设备装好之后每天怎么运行"的问题，决策变量包括各时段电解槽功率、储氢罐充放氢量、氨合成塔负荷率、并网模式下买卖电功率、离网模式下储能电池充放电功率。这些变量决定系统的运行成本和收益，属于短期决策。
容量和调度之间是强耦合的。你装多大的电解槽，直接决定了调度层能安排多少电解功率；反过来，如果调度层发现某个容量配置下运行成本特别高，上层就应该调整容量方案。所以必须把两层问题嵌套在一起求解，也就是常说的"容量-调度联合优化"。我这次复现采用的是让Cplex一次性求解整个混合整数线性规划模型，没有做外层启发式和内层LP的嵌套迭代，这样处理的好处是能保证全局最优，坏处是变量和约束规模明显变大，对建模规范性的要求更高。
1.3 为什么选Cplex而不是其他求解器
做这类中大规模MILP问题，Cplex几乎是我第一个想到的求解器。首先是求解性能确实扛得住，带几百个二进制变量的MILP模型，它通常在几分钟内就能收敛到1%以内的最优间隙；其次Matlab工具箱里直接带cplex接口，和YALMIP、CVX这种建模层的配合也成熟。
当然也有替代方案。Gurobi在纯LP和MILP上跟Cplex不相上下，而且许可证对学生更友好；如果模型是凸的，用CVX或者YALMIP内置的Sedumi、SDPT3也可以，但处理大规模混合整数问题时收敛速度明显不如商业求解器。另外国内常用的一些开源求解器如SCIP、CBC，能力也不弱，但在这个模型规模下求解时间会拉得比较长。
我的建议是：如果学校或课题组有Cplex授权，直接用Cplex；如果没有，优先试Gurobi，代码改动很小；实在没有商业求解器可用，再考虑SCIP，但要做好求解时间翻好几倍的心理准备。
2. 数学模型怎么建：目标函数与约束拆解
2.1 目标函数设计
目标函数采用全生命周期成本最小化，把投资成本通过资本回收系数折算到每年的等年值，再加上年运行成本，形成一个可比的经济性指标。
投资成本部分涵盖五项：风机、光伏、电解槽、储氢罐、氨合成塔。这里有一个关键的折算处理：每项设备单位投资成本乘以其容量，再乘以资本回收系数CRF，即
[ C_{inv} = CRF \times \sum_{i \in S} (c_i \times Cap_i) ]
其中CRF的计算公式是常见的那套基于折现率和项目寿命期的年金公式。如果项目寿命是20年、折现率取6%，CRF大约在0.0872左右。这个系数要把折现率取准了，取低了会低估投资的真实成本，取高了会过度惩罚初投资，导致优化结果偏向于少装设备。我在这里踩过坑，后面会细说。
运行成本分两块：并网模式下与电网交互的购电费用减去售电收益，离网模式下主要是弃风弃光损失对应的虚拟成本。有储氢罐以后，系统实际上有一定的跨时段能量转移能力，所以运行成本的计算必须按小时展开，不能简化成单个典型日的总量。
如果合成氨产品有出售收益，目标函数也可以改成净成本最小化，即总成本减去氨销售收入。取决于你研究的是"系统经济性最优"还是"给定氨产量下成本最低"，两者建模方式略有不同。我这里采用的是后者，即固定年氨产量，在这个前提下做容量和调度的联合优化。
2.2 关键约束条件
约束条件是这套模型的灵魂，也是复现时最容易出错的地方。我把核心约束按物理环节拆开来梳理。
第一条是各时段功率平衡约束。并网模式下每小时满足：风电出力、光伏出力、电网购电功率、储能放电功率之和等于电解槽耗电、电网售电功率、储能充电功率之和。注意离网模式下购售电项全部置零，储能的作用就被放大了。功率平衡约束必须严格满足，这是整个调度可行域的地基。
第二条是电解槽运行约束。电解槽不是任意功率都能运行，它有一个最小运行负荷率，一般碱性电解槽在20%到40%之间，PEM电解槽可以低到5%。这个约束用一条不等式就能表达：电解功率要么为0，要么位于最小和最大出力之间。这就引入了二进制变量，也是模型变成MILP的原因。如果忽略最小负荷率，优化器会让电解槽以极低的功率运行，这在物理上是不可能实现的。
第三条是储氢罐动态平衡约束。储氢罐在每个时段末尾的储氢量等于上一时段末尾的储氢量加上本时段产氢量减去本时段供氨合成塔用氢量，同时考虑电解槽的效率折算和储氢罐自身的损耗。这个约束是跨时段耦合的，它把调度问题从单时段优化变成了真正意义上的动态优化问题。
第四条是氨合成塔运行约束。氨合成塔有两个关键限制：设计产能上限和最低稳定运行负荷率。因为合成塔要求连续稳定运行，基本上一旦开机就不能停，停机成本极高，所以模型里通常会要求氨合成塔在整个时间范围内保持最低负荷以上运行，这样储氢罐就必须保证持续供氢。
2.3 非线性项的线性化处理
风光出力、电解槽效率这些物理关系天然都是非线性的，但Cplex求解MILP的强项在线性规划松弛，非线性约束要么手动线性化，要么引入分段线性近似，不能让求解器直接处理非线性项。
风电出力跟风速的关系是典型的分段函数：低于切入风速不发电，在切入风速和额定风速之间按三次方多项式上升，在额定风速和切出风速之间保持额定，超过切出风速停机。我对风速数据进行预处理，把每个时段的风电最大出力直接作为输入参数传入模型，这样就不需要在该时段内引入非线性约束，这是一个很实用的简化技巧。
电解槽的制氢效率随运行功率变化，严格来说是一条非线性曲线。我采用分段线性化的方式：把电解槽功率分成三段，每段定义一个等效产氢效率，然后引入二进制变量保证只能落在其中一段，用大M法表达分段激活逻辑。这样模型就变成了MILP，Cplex可以直接处理。
3. Matlab调用Cplex实现容量-调度联合优化
3.1 开发环境配置与Cplex接入
先把环境搭好。我用的是Matlab R2022a搭配IBM ILOG CPLEX Studio 12.10，操作系统是Windows 10。安装CPLEX时需要注意，务必勾选"Matlab Integration"组件，否则之后还要手动配置路径。
安装完成后，在Matlab里设置路径。最简单的方式是把CPLEX的Matlab接口目录添加到搜索路径，这一步我通常会写进启动脚本：
addpath('C:\Program Files\IBM\ILOG\CPLEX_Studio1210\cplex\matlab\x64_win64');配置好后验证一下是否生效：在命令行输入cplex，如果弹出CPLEX对象的帮助信息，说明接口已经打通。这一步卡住的人不在少数，常见的原因是路径写错或者Matlab版本跟CPLEX不兼容。CPLEX 12.10官方支持的Matlab版本是R2017b到R2022a，太新的Matlab版本偶尔会出现mex文件不兼容的情况，这时候要么换CPLEX版本，要么在配置里选择"使用旧版Matlab兼容模式"。
3.2 数据准备与参数预处理
模型跑得对不对，很大程度取决于输入数据准备的精细程度。
风资源数据采用某地全年的逐小时平均风速序列，光伏数据采用对应年份的逐小时太阳辐照度序列，这两个时间序列是全年8760小时的。直接把全年8760个小时全部放进优化模型，变量规模和求解时间都会失控，所以需要进行典型日选取。我采用的是K-means聚类方法，用风速、辐照度、气温三个特征把全年数据聚成4类，每类取聚类中心最邻近的一天作为典型日，并统计该类的天数占比作为权重。这样模型从8760小时压缩到96个小时（4个典型日乘以24小时），求解规模大幅下降，但代表性损失不大。
设备成本参数我用的是目前行业内比较认可的一组数据：风机单位投资7000元每千瓦，光伏4000元每千瓦，碱性电解槽5000元每千瓦，储氢罐按储氢量计约600元每标准立方米，氨合成塔按年产能分摊后大约每吨氨800至1200元。运维成本按设备投资的固定比例取2%。
电价数据方面，并网模式参考国内某省级电网的一般工商业分时电价：峰时段1.2元每千瓦时、平时段0.75元每千瓦时、谷时段0.35元每千瓦时，同时设置上网售电电价0.4元每千瓦时。分时电价的引入会直接影响调度层的决策，特别是电解槽的启停安排，谷电阶段会主动加大电解功率把氢气存起来供白天氨合成塔使用。
3.3 核心代码框架与求解流程
模型代码我按功能模块组织：数据加载模块、参数设置模块、变量定义模块、约束构建模块、求解模块、后处理模块。
目标函数的构建是核心，需要把所有变量的成本系数汇总成一个向量f，然后传给Cplex。变量的排序要谨慎，我在代码注释里详细记录每个变量在决策向量中的索引位置，避免后续构建约束矩阵时串位。实际调试中经常出现"约束对不上变量"的问题，排查起来非常痛苦，所以从一开始就要把索引管理好。
约束矩阵的构建采用逐条添加的方式。先初始化稀疏矩阵，然后按约束类型逐行填入系数。这里有一个经验性的建议：约束矩阵用稀疏格式存储，否则96个时段加上几十个变量，虽然内存不至于爆掉，但求解器预处理效率会明显下降。
求解代码的核心部分大致是这个流程：
% 创建Cplex模型对象 cplex = Cplex('WindPV_H2_NH3'); % 设置目标方向为最小化 cplex.Model.sense = 'minimize'; % 目标函数系数向量，按连续变量在前、二进制变量在后的顺序排列 cplex.Model.obj = f; % 不等式约束矩阵 Aineq：变量系数矩阵 % 注意同时设置 lhs 和 rhs，实现上下界约束 cplex.Model.A = A_ineq; cplex.Model.lhs = lhs_ineq; cplex.Model.rhs = rhs_ineq; % 等式约束 cplex.Model.A = [cplex.Model.A; A_eq]; cplex.Model.lhs = [cplex.Model.lhs; b_eq]; cplex.Model.rhs = [cplex.Model.rhs; b_eq]; % 变量边界 cplex.Model.lb = lb; cplex.Model.ub = ub; % 变量类型：C代表连续，B代表二进制 cplex.Model.ctype = ctype; % 求解器参数设置 cplex.Param.timelimit.Cur = 1800; cplex.Param.mip.tolerances.mipgap.Cur = 0.005; cplex.Param.mip.strategy.startalgorithm.Cur = 4; cplex.Param.threads.Cur = 8; % 执行求解 cplex.solve(); % 获取结果 if cplex.Solution.status == 101 || cplex.Solution.status == 102 x_opt = cplex.Solution.x; obj_opt = cplex.Solution.objval; else error('求解未成功，状态码：%d', cplex.Solution.status); end状态码101和102分别对应"最优解"和"达到最优间隙的可行解"，这两个状态都可以接受。求解完以后，把决策变量按之前的索引规则还原成各个设备的容量和逐小时功率数据，再跑一个校验脚本，把所有约束重新检查一遍，确认无违规量。
3.4 求解性能调优的实用策略
Cplex求解MILP的速度跟模型质量高度相关，同样的数学问题，建模方式不同求解时间可能差出一个量级。
我调试过程中最见效的一个优化手段是给二进制变量设置合理的初始解。提前跑一遍启发式策略，把容量变量固定到几个合理的经验值上，然后只优化调度层，得到一个初步调度方案，再把这个方案中的二进制变量取值作为MIP start传给Cplex。这个操作可以显著削减分支定界树的搜索空间。
另一种策略是适度放松MIP间隙容忍度。学术研究追求1e-4的间隙当然好，但工程实践中，0.5%到1%的最优间隙已经足够支撑决策分析。在代码里设mipgap为0.005，求解时间通常会从几十分钟缩短到几分钟，而目标函数值的差异往往不到0.3%。
还要注意大M参数的取值。处理"电解槽要么关要么在最小功率以上运行"这类逻辑约束时，大M的取值不能太大，大M太大容易导致数值病态，求解器在判可行性时会出问题；大M太小又有可能把可行域错误地截断。我的经验是M取"该设备全年最大运行功率对应的影响量级乘以1.5"就够了。
4. 并网 vs 离网：容量配置结果与调度规律
4.1 典型日的选取与验证
典型日聚类完成后，先做一步验证工作：把聚类得到的四个典型日按权重叠加，算出来的全年总发电量和全年总制氢量跟直接用8760小时逐时数据计算的结果做对比。偏差在5%以内，说明这四个典型日确实能代表全年的资源特性。
我复现时最终选出的四个典型日大致对应冬季大风日、春季过渡日、夏季高温高辐照日、秋季多云日。冬季大风日风速高但辐照弱，光伏出力小；夏季光伏出力大但风速相对较低。这种差异化组合正好能覆盖风光的互补特性，也让优化器有足够的场景去权衡容量配比。
4.2 两种模式下的最优容量配置对比
求解完成后，最直观的结果是并网和离网两种模式下系统最优容量配置的差异。为了对比统一，我在两种模式下都设定相同的年氨产量目标。
表：并网与离网模式下最优容量配置对比
| 设备 | 并网模式 | 离网模式 |
|---|---|---|
| 风电装机容量（兆瓦） | 58 | 76 |
| 光伏装机容量（兆瓦） | 42 | 66 |
| 电解槽额定功率（兆瓦） | 51 | 62 |
| 储氢罐有效容积（万标准立方米） | 12.8 | 21.5 |
| 储能电池容量（兆瓦时） | 无 | 38 |
并网模式下系统大幅依赖电网购电来补充短时功率缺额，风电和光伏的装机容量相对较小，储氢罐也不需要做得很大，因为电网本身就是一个无限的缓冲池。离网模式下风电光伏装机容量分别提高了31%和57%，储氢罐扩容68%，这是因为离网系统必须靠源端冗余和储氢来应对风光出力的连续波动，没有任何外部支撑。
从经济性指标来看，离网模式的单位氨生产成本比并网模式高约15%到20%。多出来的成本主要来自更大的初投资，以及为了满足"零外购电"约束而必然出现的弃风弃光电量。这一结果符合预期，但量化出来还是很有参考价值的：如果项目选址附近有稳定电网接入条件，优先选并网模式，单位成本明显更低；只有定位在零碳示范项目或者偏远无电网地区，才值得考虑离网方案。
4.3 调度策略与经济性分析
调度结果的规律比容量结果更有意思。并网模式下，电解槽的运行曲线明显跟着电价走：谷电时段电解功率爬升到接近额定功率，峰电时段则降功率甚至停机，氢气需求缺口依靠储氢罐中预先存好的氢气来补足。这实际上就是利用氢储能的跨时段转移特性做电价套利。
离网模式下没有电价信号可以跟随，调度策略完全由风光出力和储氢罐状态决定。典型大风日的调度曲线显示：风电出力高峰时段电解槽满负荷运行，同时储能电池充电；风电出力低谷时段储能电池放电维持电解槽基本负荷，氨合成塔始终保持不低于70%的设计负荷运行。储氢罐在这个场景里扮演的角色就是能量缓冲池，风大的时候囤氢，风小的时候放氢，保证合成氨连续生产。
经济性分析方面，把成本拆开看：并网模式下运行成本中购电费用占比接近40%，这是降本的关键突破口；离网模式下投资成本占比超过75%，运行成本主要是设备折旧摊销，边际运行成本很低，多产一吨氨的增量成本远低于并网模式。所以长期来看，如果风光设备成本持续下降，离网模式的经济性劣势会逐步收窄。
5. 调参、踩坑与排查实录
5.1 求解失败与可行性问题
复现过程中最难受的问题不是求解慢，而是模型不可行。Cplex返回infeasible的时候，第一反应不要怀疑求解器，先怀疑自己的约束建模。
我遇到的一次典型问题是：氨合成塔的最低运行负荷约束和储氢罐的初始储氢量约束冲突了。模型假设第一天8点开始运行，储氢罐初始储氢量设得偏低，而氨合成塔要求全天候最低负荷运行，导致晚上氢气供应不足，约束冲突。排查方法是用Cplex的冲突定理解析功能，让求解器自动找出一组不可行约束子集，很快就能定位到问题所在。
cplex.conflict.minConflict(); cplex.conflict.get();这个功能在调试初期特别实用，但需要提醒的是，对大规模模型启用冲突检测会显著增加求解时间，建议单独跑一个调试版的模型文件来做冲突分析，不要在生产版本里开这个功能。
5.2 数据量纲和精度问题的几个典型坑
量纲问题是新手最容易踩的坑，我也栽过一次。风速数据用的单位是米每秒，功率用兆瓦，储氢量用公斤，而电价用的单位是元每千瓦时，几个量纲混在一起，约束矩阵的系数跨度从1e-3到1e5，数值范围差异过大，Cplex内部处理时容易出现数值警告。
解决方案是统一量纲体系：所有能量统一用兆瓦时，氢量统一用吨，成本统一用万元。功率平衡约束里所有项都变成兆瓦，物料平衡约束里所有项都变成吨，这之后系数范围就基本控制在1e-2到1e3之间了，数值上的报警信息也消失了。
另一个坑是储氢罐动态方程里的时间步长。我们用的都是1小时步长，但电解槽制氢量的计算涉及连续时间到离散时间的转换。用1小时步长的时候，某时段电解槽以50兆瓦功率运行，对应的氢气产量是该时段总制氢量，这没有歧义。但如果后期把模型改成15分钟调度，就有必要核查效率系数是否需要按时间步长重新换算，稍不留神就会导致全年的氢气平衡不闭合。
5.3 性能瓶颈与后续扩展建议
模型规模增大后，求解时间会非线性上升。96个时段的基础模型，Cplex大约在5分钟内跑到0.5%间隙；如果把典型日扩展到12个，求解时间会涨到40分钟级别。如果后续想把8760小时全部纳入模型做精细化分析，建议采用Benders分解或者拉格朗日松弛把原问题拆成容量主问题和调度子问题，用迭代的方式收敛，而不是让Cplex硬扛。
从模型扩展角度看，几个方向值得考虑。第一，把电解槽的动态响应特性加入约束，即电解功率的爬坡速率限制，这会增加约束数量但更贴近实际；第二，考虑氨合成塔的启停成本，在目标函数里引入固定启动费用，这样调度优化会在"维持低负荷运行"和"停机再启动"之间做更真实的权衡；第三，加入可再生能源出力不确定性的随机优化或鲁棒优化扩展，这是目前学术研究的热点方向，也是工程实际中必须面对的问题。
我这次复现的全部代码和数据文件按照模块化方式组织，每个模块都有详细的注释，跑完一遍以后，最大的感受是：这类问题的难点从来不在求解器本身，而在于把物理过程转化为数学约束时做的那些简化和取舍。每一个线性化近似、每一个典型日聚类、每一处效率系数的取值，都会影响最终容量配置的经济性和可行性。建模的时候多想一步，调试的时候就能省十步。
