前面一直有朋友私信我问综合能源方向的学生项目到底该怎么落地。说实话这几年只要做过综合能源系统建模几乎都绕不开同一个问题在碳约束和可再生能源高比例接入的双重压力下传统的热电联产CHP调度模型已经不够用了必须把电转气P2G和碳捕集系统CCS也纳入同一个优化框架里。但这个纳入说起来容易做起来远比想象中复杂——电、热、气、碳四条平衡线互相牵扯光靠直觉或者静态计算根本算不清最后我还是回到Matlab里老老实实建模、做24小时滚动优化才把这个问题彻底想明白。这篇内容我准备讲透三层东西一是P2G和CCS接入CHP系统后设备模型分别长什么样、参数怎么定二是整个优化问题的目标函数和约束条件如何在Yalmip里写出来、用Gurobi求解三是我在调试过程中踩过的几个大坑包括量纲混乱、非线性项没法收敛、热负荷约束导致无解这些新手必踩的问题。适合正在做低碳园区、微能源网、综合能源优化调度的研究生或工程技术人员作为参考代码思路和参数直接可以迁移到自己的项目里。1. 为什么非要把P2G和CCS绑在CHP上这不是设备堆叠而是碳回路的闭环1.1 单纯CHP调度模型面临的碳账单先说传统CHP的调度。无论是背压式还是抽凝式传统模型的逻辑都很简单给定电负荷和热负荷曲线在热电耦合可行域内找一个最低燃料成本的工作点。目标函数里只有购气费用和购电费用约束只有功率平衡、热功率平衡和设备容量。但一旦引入碳排放交易机制事情就变了。排碳变成真金白银的成本CHP机组每发一度电、每供一单位热都得付相应的碳价。这时候调度模型里就出现了一个新的自由度——机组是维持高负荷运转、买了碳配额换取更多电热产出还是适当降负荷、少排碳但去现货市场买电这笔账不是拍拍脑袋能算出来的必须放进优化目标里一起算。更关键的是CHP烟气里的CO2浓度比普通燃气锅炉还要高这恰恰是碳捕集最经济的场景。于是CCS自然就成了下一个要加进来的子系统。1.2 P2G给捕获的CO2找了一条出路原本加了CCS之后捕获的CO2要么压缩后地质封存要么用于驱油但在园区级综合能源系统里这两种去向都不太现实——没有封存场地也找不到买家。这时候P2G提供了一个非常漂亮的解决方案用可再生能源富余电力电解水制氢再用氢气和CO2通过甲烷化反应合成天然气SNG。注意看这个链条CHP排碳 → CCS捕集 → P2G甲烷化消耗碳 → 产出的SNG又重新回到燃气系统里供给CHP。碳不往系统外排而是在内部循环。这个碳回路是整套系统最精彩的地方但也是建模最麻烦的地方。麻烦在两点第一甲烷化反应对CO2的消耗量与CCS的捕集量在同一时刻需要平衡如果捕多了没地方存捕少了P2G就得停机第二CCS的再生能耗、P2G的耗电量都会反过来拉低CHP的净电出力形成一个复杂的能-碳双向耦合。1.3 模型最终要回答的四个实际问题其实一个课题值不值得做就看它能不能回答几个原先回答不了的问题。这套P2GCCSCHP模型我认为核心要回答四个在给定的风电、光伏、电负荷、热负荷曲线以及分时电价和碳价条件下CHP电出力、热出力、CCS捕集率、P2G功率、储气/储热充放分别取多少系统总运行成本最低CCS捕集率怎么跟着碳价走碳价高的时候多捕电价高的时候要不要为了省电而放弃捕集储气罐和储热罐的容量配置对系统成本影响有多大最优24小时充放策略是什么加入P2G之后弃风率降了多少天然气净购入量是升了还是降了这四个问题就是整个Matlab程序的设计目标。我后面每一个章节本质上都是在为回答这些问题做准备。2. 先画系统拓扑再写代码能流、碳流、热流的完整连接关系2.1 整套系统的物理拓扑我自己写代码之前习惯先画一张图把设备接口关系理清楚这比直接开写约束条件高效得多。这套系统的物理结构如下能源输入侧外部电网从市场购电、天然气网购气、风电场和光伏阵列。转换设备CHP热电联产机组、电解槽、甲烷化反应器、碳捕集装置、电锅炉作为储热的补充热源。存储设备储气罐储存SNG、储热罐储存热水的显热。负荷侧电负荷、热负荷。能量路径要分三条线看。电的路径是外购电风电光伏CHP发电 → 电负荷CCS耗电P2G电解槽耗电电锅炉耗电。热的路径是CHP热出力电锅炉热出力 → 热负荷储热罐蓄热CCS的再生热如果需要热也从这条线取。气的路径是天然气网购气P2G产SNG → CHP燃料消耗储气罐充放。碳流的路径则是CHP烟气CO2 → CCS吸收塔 → 解吸 → CO2储罐 → 甲烷化反应器 → SNG → 回到天然气侧。这是一条真实的物质循环链在建模时需要注意CO2储罐的SOC与储气罐SOC是两个不同的状态变量不可以混在一起。2.2 四条平衡线要同时闭合整个优化模型的核心就是让以下四条平衡约束在每个时段都同时成立电功率平衡P_buy P_wt P_chp_e - P_ccs_e - P_p2g - P_eb P_load热功率平衡H_chp H_eb H_dis_tank H_load H_charge_tank H_ccs气平衡G_buy G_sng G_chp G_charge - G_dis碳平衡E_chp_total E_capture E_net净排入大气的碳量 E_p2g_consume进入甲烷化的碳量注意碳平衡方程是最容易出错的地方。E_chp_total是完全由燃料燃烧产生的CO2总量E_capture是CCS捕集掉的量E_p2g_consume是P2G甲烷化反应消耗掉的CO2。捕集量和消耗量之间不是恒等关系如果E_capture E_p2g_consume多出的CO2有可能进入储罐储存如果E_capture E_p2g_consume就需要从外部补充碳源比如购买工业CO2这一点在实际工程里常被忽略。碳交易市场计算的净排放量用的是E_chp_total - E_capture这个净口径即实际排入大气的碳而不是燃料燃烧产生的总碳量。我自己一开始就在这里犯了错把总排碳直接算成碳排放结果碳交易成本多了将近一倍折腾了半天才找到原因。2.3 为什么要设置三种对比运行模式为了把每个设备的价值拆开看代码里我实现了三个模式M1纯CHP无CCS、无P2G、无储气。这是基准模式作为成本下限参考。M2CHPCCS捕集CO2但不接入P2G碳捕集后的CO2视为外运或封存。M3CHPCCSP2G储气储热全配置碳回路完全闭环。三种模式跑同一套负荷与电价数据对比结果就能厘清CCS单独带来的成本压力有多大P2G补上之后又省了多少弃风损失。这也是写论文时最有说服力的一组对比图。3. 核心设备建模给每个模块一件数学外衣3.1 CHP机组的热电耦合可行域CHP建模的核心不是效率公式而是热电可行域。背压式CHP的电出力完全由热出力决定P c * H灵活性几乎为零。而抽凝式CHP的热电可行域是一个凸多边形由几条线性不等式围出来0 ≤ P_chp ≤ P_max 0 ≤ H_chp ≤ H_max P_chp ≥ k1 * H_chp b1 P_chp ≤ k2 * H_chp b2k1、k2分别是可行域上下边界的斜率。我取的典型参数是P_max 200 MW、H_max 160 MW、k1 0.15、k2 0.85这代表机组在最大凝汽工况和最大抽汽工况之间的调节范围。也就是说在同一热出力下机组电出力不是唯一的而是可以在一个区间内调节——这正是优化调度能够发挥作用的自由度。CHP的燃料消耗严格来说与电、热出力之间是二次非线性关系。但在MILP混合整数线性规划框架下二次项没法直接处理工程上最常用的做法是分段线性化。把热出力轴切成8到10段每一段内燃料消耗按线性函数近似。分段数太少可行域顶点附近的工况点误差会特别大分段数太多求解时间又会显著拉长。实测下来10段是一个性价比比较合适的取值。3.2 碳捕集系统能耗模型和捕集率边界CCS我用的是燃烧后胺法捕集的简化稳态模型。捕集对象是CHP烟气关键变量是捕集率α0到0.9之间的连续变量也可离散为0/1启停。捕集量Q_cap α × E_co2_prod其中E_co2_prod是CHP燃料燃烧产生的总CO2量按每MWh天然气燃料约产生0.2吨CO2估算。CCS的能耗要分两笔账电耗约每吨CO2耗电0.060.10 MWh用于压缩机和泵组E_ccs_e β × Q_cap。热耗约每吨CO2耗热0.20.4 MWh用于胺液再生H_ccs_th γ × Q_cap。这两笔能耗必须在电平衡和热平衡中同时体现很多人只算了电耗却漏了热耗结果导致储热罐容量设计偏小。热耗是我在代码里特意单独拆出来的一条约束变量这个细节帮我在后面做敏感性分析时发现了一个有趣的现象碳价上升时CCS捕集率上升但CCS的热耗也随之增加导致CHP原本供给热负荷的热量被挤占储热罐的放热深度明显加大。如果没有储热缓冲系统的热供应很快就会出现缺口。捕集率还有个上限问题。实际胺法系统吸收塔负荷不可能无限提升通常上限取0.850.90超过这个数值后单位能耗会急剧增加。我在模型里直接用硬上限0.9处理如果要做更精细的分析可以把能耗系数改成捕集率的递增函数但这会引入非线性项后面会讲到对应的线性化处理办法。3.3 电转气系统电解槽和甲烷化单元的建模关键P2G分两级建模。第一级是电解槽效率η_el通常在0.60.75之间。考虑实际运行约束电解槽不允许极低负荷运行我设了20%的最小运行负荷低于这个值必须停机或者完全关闭。这一条需要二进制变量如果不用二进制变量优化器会给出在5%负荷运行电解槽这种物理上不成立的结果。第二级是甲烷化反应器。反应式为4H2 CO2 → CH4 2H2O。从质量守恒的角度每生产1 MWh的SNG需要消耗约0.2吨CO2。甲烷化过程本身是放热反应温度通常在250400摄氏度从动态响应的角度来说它不宜频繁启停我给它加了一个最小连续运行时间的约束。这同样需要二进制变量是模型从纯LP变成MILP的关键原因。补充一点P2G设备在系统中的地位其实是电力弹性负荷。Wind和光伏出力大的时候弃风弃光已经发生或预测即将发生时P2G以低价甚至负价的地区消纳这部分电力。但作为代价P2G运行起来后系统电负荷会突然多出一大块这个冲击要靠外部电网和CHP的调节能力去吸收。如果P2G容量设计过大夜间风电高峰时段的电平衡压力会很明显这就是为什么必须做优化而不是简单能开就开。3.4 储气罐与储热罐的动态约束储气罐和储热罐属于同一个数学框架——能量存储方程SOC(t1) SOC(t) P_in(t) - P_out(t)其中SOC是罐内能量存量的归一化值P_in、P_out是充放功率。边界约束包括容量约束0 ≤ SOC(t) ≤ C_tank充放速率约束P_in(t) ≤ R_in_maxP_out(t) ≤ R_out_max周期约束SOC(24) SOC(0)保证一个调度周期内储能量回归初始值便于滚动调度衔接储热罐的价值尤其值得强调。CHP加入CCS后热耗变大热出力灵活性下降储热罐能够把CHP的热出力与热负荷在时间上解耦让机组在电价低谷时降低电出力同时热出力也降但依靠储热满足热负荷电价高峰时追求满发、把过剩热存起来。热电解耦是这套系统经济性提升的主要来源。如果项目里没有储热P2G和CCS的加入只会让CHP更加僵化优化结果几乎没有任何意义。4. 优化模型与Matlab实现YalmipGurobi从零搭建4.1 目标函数把成本账列清楚目标函数是系统24小时总运行成本最小化单位为元。五个分项购气费用G_buy沿时间累加乘天然气价格元/MWh购电费用P_buy乘分时电价元/MWh峰谷价差足够大的情境下这部分是优化空间最大的项碳交易费用(E_total_net - E_quota) × 碳价元/吨E_net为正则买碳配额为负则可出售盈余配额获得收益设备运维费用按出力线性比例估算弃风惩罚为了让模型尽量消纳风电对弃风量给定较高的惩罚单价元/MWh需要注意的是碳价和弃风惩罚这两个参数对结果影响极其敏感。如果弃风惩罚设得太高优化器可能会牺牲经济性强行消纳风电如果设得太低P2G根本不会被启用碳回路也就形同虚设。4.2 约束条件的Yalmip等价写法我用的建模框架是Matlab Yalmip Gurobi。决策变量定义如下P_chp sdpvar(1, 24); % CHP电出力MW H_chp sdpvar(1, 24); % CHP热出力MW alpha sdpvar(1, 24); % CCS捕集率0~0.9 P_p2g sdpvar(1, 24); % P2G耗电MW SOC_h sdpvar(1, 24); % 储热罐能量状态MWh SOC_g sdpvar(1, 24); % 储气罐能量状态MWh u_chp binvar(1, 24); % CHP启停状态 u_p2g binvar(1, 24); % P2G启停状态约束组装的核心思路是逐时段写约束C []; for t 1:24 % 热平衡CHP热出力 电锅炉 储热放热 热负荷 CCS热耗 储热吸热 C [C, H_chp(t) H_eb(t) SOC_h(t)/eta_h_dis ... H_load(t) H_ccs(t) SOC_h(t-1)*eta_h_chg H_eb_loss(t)]; % 电平衡外购电 风电 CHP电出力 电负荷 CCS电耗 P2G耗电 电锅炉耗电 C [C, P_buy(t) P_wt(t) P_chp(t) ... P_load(t) P_ccs_e(t) P_p2g(t) P_eb(t)]; % CHP热电可行域 C [C, P_chp(t) k1 * H_chp(t) b1]; C [C, P_chp(t) k2 * H_chp(t) b2]; % CCS捕集量与能耗 C [C, Q_cap(t) alpha(t) * (co2_rate * F_chp(t))]; C [C, P_ccs_e(t) beta * Q_cap(t)]; C [C, H_ccs(t) gamma * Q_cap(t)]; C [C, alpha(t) 0.9]; % P2G产气量与CO2消耗 C [C, P_sng(t) eta_p2g * P_p2g(t)]; C [C, CO2_p2g(t) co2_per_sng * P_sng(t)]; % 碳平衡净排碳 总排碳 - 捕集碳 - P2G耗碳 C [C, E_net(t) co2_rate * F_chp(t) - Q_cap(t) - CO2_p2g(t)]; end这段代码是高度浓缩的核心逻辑。实际运行中还需要配合储气罐的SOC递推、CO2储罐的SOC递推、储热/储气的容量约束这些用同样的方法累加进C矩阵即可。注意SOC_h(0)需要预先给定初值否则Yalmip会报未定义变量错误。4.3 求解器选型为什么是Gurobi而不是自带的linprog一旦引入了0-1变量启停、分段线性化问题就变成MILP。Matlab自带的linprog只能处理纯线性LP是没办法直接求解MILP的。我选择的组合是Yalmip做建模层、Gurobi做求解层。Gurobi对MILP的处理能力显著优于开源求解器尤其是大规模整数规划问题差距很大。求解器配置有几个关键参数要调ops sdpsettings(solver, gurobi, ... gurobi.MIPGap, 0.001, ... gurobi.TimeLimit, 300, ... verbose, 1); result optimize(C, objective, ops);MIPGap整数规划相对间隙设0.1%已经足够工程精度设为0会导致极端情况下求解时间急剧膨胀。TimeLimit设300秒是防止某个工况点长期卡在分支定界树里。我遇到过Gurobi在200个整数变量时跑了一个多小时的情况就是因为MIPGap设成了0。4.4 结果输出与绘图一张图讲清楚一套数据优化结束后我习惯用yyaxis双轴绘制电功率平衡曲线和热功率平衡曲线用堆叠面积图展示电负荷中外购电、风电、CHP发电各自的贡献再用一张单独的CO2曲线子图展示总排碳、捕集、P2G消耗、净排放四条曲线。figure; subplot(2,2,1); stairs(1:24, P_wt, k-); hold on; stairs(1:24, P_chp, r-); ...所有结果同时导出为.xlsx方便后续做敏感性分析时反复取数。曲线绘制的一个细节提醒不同变量的单位在数值上可能巧合相等比如1 MW × 1h 1 MWh但画图时Y轴一定标注清楚否则汇报时候很容易被问住。5. 典型日案例与结果解读三种模式的曲线说明了什么5.1 案例数据的设定逻辑为了让结果有说服力我构造了一个北方冬季典型日案例热负荷白天高、夜间略降风电出力夜间大、白天小电负荷呈双峰曲线。分时电价采用峰平谷三段式峰时电价是谷时的3倍。碳价设为100元/吨。CHP机组额定容量200 MW配30 MW电锅炉P2G容量20 MWCCS最大捕集能力对应0.9捕集率储热罐容量80 MWh。5.2 M1基准模式热跟着负荷走电出力很僵硬M1模式下CHP没有储热、没有CCS机组只能紧紧跟随热负荷运行。夜间热负荷虽然不算特别高但为了保证供热CHP电出力被强制推高而此时风电大发、电负荷又低多余的电只能弃掉。这个场景特别典型冬季供热季的热电矛盾本质上是强迫出力造成的风电消纳困难。M1的弃风率大概在18%22%之间这个数字很容易通过模型复现。5.3 M3联合模式夜间风电转化为气并重新参与供给M3模式的结果明显不同。夜间风电大发时段优化器给出的策略是提高P2G功率用电解槽吃下富余风电同步提高CCS捕集率把CHP烟气中的CO2尽量捕下来喂给甲烷化反应器CHP热出力维持较高水平多余热量存进储热罐。于是CO2平衡图上可以看到三条曲线在夜间出现了明显的汇聚总排碳、捕集碳、P2G耗碳三条线在夜间同时达到高峰而净排碳曲线反而下探。更直观的是储热罐SOC曲线白天电价高峰时段CHP满发与高捕集率并行储热罐放热夜间电价低谷时段储热罐蓄热CHP在保热电出力下运行。储热罐的削峰填谷作用一目了然。M3模式弃风率降到5%以下天然气购气量则因SNG的补充而下降了约15%。5.4 碳价的敏感性拐点不是越低碳价越好我额外做了一组碳价从50元/吨到500元/吨的敏感性扫描。结果有个明显的拐点碳价从50涨到200元/吨时CCS捕集率和P2G出力几乎线性上升净碳排显著下降但碳价200元/吨之后再提高碳价总成本的降幅趋缓甚至出现成本反弹。原因是CCS能耗和P2G能耗带来的购电成本增量与省下的碳成本逐渐打平。换句话说这套系统的碳减排存在一个经济最优区间碳价一味提高未必让系统总成本继续下降。这个结论对政策制定者是有参考价值的对做工程的朋友来说也意味着设备容量和运行策略要跟着碳价水平做适配而不是一次性设定死了。6. 实操过程踩过的坑与调试经验6.1 量纲混乱CO2用吨还是用MW差44倍的教训我最早一版代码里CCS模型用吨CO2P2G用kgCO2结果耦合约束E_net co2_rate × F_chp - Q_cap - CO2_p2g怎么都对不上。检查了整整一天最后发现CO2_p2g算出来比Q_cap大两个数量级才意识到一个用了吨、一个用了kg。这个问题的隐蔽性在于如果只是倍数关系模型依然能求出一个最优解只不过物理意义完全错误。实操建议是全程统一能量和质量单位并且在每个变量定义处用注释标清单位。6.2 非线性项导致的求解失败甲烷化反应的CO2消耗量原本是温度的非线性函数我第一版想用更精确的模型结果Yalmip报非凸二次约束无法求解Gurobi直接拒绝跑。后来改用线性关系CO2_p2g co2_per_sng × P_sng把反应效率当作常数问题立刻解决。对工程级调度模型来说在线性框架内保留核心物理规律比一味追求精度但模型不可解要重要得多。6.3 热负荷约束太紧的infeasible问题加入储热罐之后有一段时间模型频繁报infeasible。排查发现是热平衡方程写得太紧H_chp H_dis H_load H_ccs H_chg但夜间热负荷低谷时段CHP最小技术出力对应的热出力已经超过这个值导致约束无法满足、系统无解。解决办法有两个一是加入电锅炉作为热松弛二是给热平衡加一个小松弛变量并惩罚。工程上这两种做法都合理也能帮助判断系统容量配置是不是有短板。6.4 给新手的代码组织建议从1小时模型开始最后分享一个我自己反复验证有效的调试策略先跑1小时的单周期版本把所有约束、变量、参数都跑通再扩展到24小时。单周期模型里储罐的SOC递推可以简化为一个代数方程排查逻辑错误的速度快很多。之后每加入一个新设备模块先单独验证该模块的约束是否满足再接入系统级耦合。不要一上来就写800行完整代码一旦出问题排查成本会非常痛苦。说到经验层面最大的体会是模型的可解释性比精度更重要。Gurobi算出来的解再漂亮也要先回到物理本质上检查夜间弃风时段P2G有没有启动、储热罐的充放方向是否符合电价逻辑、CO2三条曲线有没有互相矛盾。如果某些环节违反直觉多半是约束写错了而不是优化器算错了。这个检查习惯能帮你省下大量排错时间也是从会跑代码走向真懂系统的必经一步。