虚拟电厂优化调度不是新话题，但把电转气、碳捕集和垃圾焚烧放在一起之后，模型的复杂度完全不一样。这段时间我在MATLAB里复现了一套这样的日前优化调度程序，用YALMIP建模、调用MILP求解器求解，今天把建模思路、代码结构和仿真结果整理出来，给同样在搞虚拟电厂调度或者减排路径研究的朋友做个参考。
这套程序的核心问题可以概括成一句话：一个虚拟电厂运营商，手里有日处理能力500吨的垃圾焚烧发电厂、40 MW风电场、20 MW光伏电站，还考虑加装一套烟气碳捕集系统和一套电转气装置。从今天零点到明天零点，该怎么安排各路电源出力、碳捕集启停、电转气产气节奏，才能在满足电力负荷和气体交付合同的前提下，让运行成本和碳排放总量都尽可能低？听起来像常规经济调度，但加入碳捕集和电转气之后，电、碳、气三条流要同时平衡，很多在传统机组组合模型里不会出现的耦合约束就冒出来了。
1. 为什么要把电转气、碳捕集和垃圾焚烧放进同一个虚拟电厂
1.1 垃圾焚烧发电在调度中的特殊地位
垃圾焚烧发电厂平时在调度模型里常被当成“可调基荷”处理，这个理解只说对了一半。垃圾进厂之后不能长期堆存，焚烧线要连续运行，向上游城管和环卫部门承诺的垃圾消纳量是硬约束，所以焚烧量的波动幅度其实很小。但在汽机侧，通过调整抽汽量，发电功率可以在一定范围内变化，并非完全刚性。
从减排角度看，垃圾焚烧厂又是碳排放大户。按常见垃圾组分估算，每吨生活垃圾焚烧产生的化石源CO2大约有0.11～0.14吨。日处理500吨的厂，每天化石源CO2排放就在55～70吨左右，全年接近2万吨。如果只看电力调度，它会稳定地发基荷电；但一旦接入碳交易机制，它的排碳成本会成为整个虚拟电厂运营中不可忽略的一块。
我复现时把垃圾焚烧出力做了简化处理，把它和入炉垃圾量线性挂钩：焚烧量决定主蒸汽产出，汽机再按热电比把蒸汽分配成功率和抽汽。这样既保留了它的“垃圾处理刚性”，又让汽机侧留出调度灵活性，模型不至于因为非线性太强而难以求解。
1.2 电转气与碳捕集之间天然存在耦合
电转气的常见路径是电解水制氢，再把氢和CO2甲烷化，得到可并入天然气管网的合成天然气。问题在于甲烷化需要稳定供应CO2，而从外部购买高纯度CO2成本很高。碳捕集系统刚好能从垃圾焚烧烟气里分离出浓度较高的CO2，两条技术路线拼在一起，恰好形成“捕集—消纳”闭环。
但在调度层面，这个闭环有一个核心矛盾：风电、光伏高峰时段电价低，适合给电转气和碳捕集系统供电，而这两套系统的功率需求会吃掉一大块风光出力；可甲烷化装置在低负荷时段刚刚启动时，碳捕集系统供给的CO2流量又未必跟得上。反过来，如果碳捕集全开而电转气停机，捕集下来的CO2没有下游去向，要么加压液化外运，要么临时储存，都会增加成本。
所以，电转气和碳捕集必须联合调度，不能分开优化。这也是我复现这个模型时觉得最有价值的地方：它会逼着你去把CCS的烟气约束、CO2储罐容量、P2G的启停约束、甲烷销量约束放在同一个优化问题里。
1.3 复现程序要解决的输入数据准备
这套程序在数据组织上分成六个模块：负荷曲线、风电功率、光伏功率、垃圾入炉量、设备参数、价格系数。风电和光伏曲线我按典型日场景构造，负荷采用工业园区典型日数据，垃圾入炉量设定为每小时20～24吨之间波动，总处理量控制在500吨/日。
设备参数按工程常见量级标定：风电40 MW、光伏20 MW、垃圾焚烧装机15 MW，碳捕集捕获率90%，捕集能耗0.25 MWh/tCO2；电转气额定功率5 MW，电解效率62%，甲烷化效率85%。所有参数不一定对应到某个具体工程，但数量级都非常接近实际设备，用来研究调度规律是足够的。如果要用到真实项目，只需要替换数据文件里的曲线和参数。
2. 复现调度模型的目标函数与运行约束
2.1 目标函数：成本、碳收益与惩罚项怎么摆平
优化目标我写成24小时总运行成本最小，拆成五个部分：
- 垃圾燃料与辅助材料成本，包括活性炭、尿素、脱硫剂等耗材；
- 各设备运行维护成本，按出力比例折算；
- 碳交易净成本，考虑免费配额和CCS捕集量抵减；
- 购气成本减去售气收入，P2G产出的甲烷按天然气价格出售；
- 弃风弃光和失负荷惩罚项，用松弛变量实现。
公式写成：
min Σt [ C_fuel(t) + C_om(t) + C_carbon(t) + C_gas(t) + C_penalty(t) ]
碳交易这一项，我的处理方法是：
C_carbon(t) = λCO2 × ( E_emit(t) − E_quota(t) − E_cap(t) )
其中E_emit是焚烧产生的化石源CO2排放，E_quota是免费配额，E_cap是CCS捕集量。λCO2作为碳价参数，设为60元/吨。如果捕集量超过配额缺口，这个值可以为负，代表碳配额结余出售收益。这里做了简化，实际碳市场中跨期存储配额会有限制，但复现调度规律够用。
惩罚项我单独拎出来说。弃风惩罚设为50元/MWh，失负荷惩罚设为2000元/MWh。前者用来鼓励系统尽量消纳新能源，但允许在极端时段少量弃风；后者是硬性成本，程序不会为了省燃料费而牺牲供电可靠性。
2.2 关键设备约束的建模细节
电力平衡约束是整条模型的主干：
P_wind(t) + P_pv(t) + P_wte(t) + P_CCS(t) + P_P2G(t) + P_shed(t) = P_load(t) + P_dump(t)
注意P_CCS和P_P2G的方向：它们本质上是负荷，所以在平衡方程中我用括号把“可用发电减去内部用电”写成净出力：
P_net(t) = P_wte(t) + P_wind(t) + P_pv(t) − P_ccs_load(t) − P_p2g_load(t)
再加上买卖电量和失负荷松弛量，等于负荷。这样写的好处是，程序调试时能清楚地看到每个时段CCS和P2G消耗了多少电。
垃圾焚烧的约束分三层：一是入炉垃圾总量的日约束，保证500吨全部处理完；二是焚烧功率上下限和爬坡约束；三是汽机发电功率与垃圾焚烧量的线性关系。碳捕集系统则增加了一层“CO2来源约束”，也就是捕集量不能超过焚烧烟气中的CO2流量乘以捕获率，同时，捕集系统的启停用二进制变量控制，冷启动状态不允许瞬时满载运行。
电转气分电解、甲烷化、储气罐三个环节建模。电解功率有上下限和爬坡约束，产氢量等于电解功率乘以电解效率；甲烷化消耗氢气和CO2，产气速度受甲烷化反应器容量限制；储气罐的容量状态变量要满足递推关系：
S_gas(t+1) = S_gas(t) + Q_produce(t) − Q_sell(t)
同时设置首末状态相等，防止模型通过“凭空消耗一天前的储气”来取巧求解。
2.3 非线性项的线性化处理
这个模型如果直接写成非线性规划，YALMIP能建模，Gurobi却没法直接解。我把所有非线性项都做了线性化处理。
碳交易部分采用分段线性碳价。现实中碳价超过一定阈值后，边际减排成本上升，碳配额交易价格往往看涨，我把它处理成阶梯状线性函数，引入两个连续变量和一个二进制变量，用Big-M法实现分段切换。
风电光伏的随机性用确定性等价场景替代。也就是先用场景法生成24个可能的风光出力场景，再用快速前向选择法缩减成3个典型场景，模型变成多场景两阶段优化：第一阶段决定机组启停和CCS启停，第二阶段在各个场景下调整功率。为了保证求解规模可控，我把启动变量设为一阶段变量，连续功率变量设为二阶段变量，这样总变量数大约800个，二进制变量约200个，MILP求解压力不大。
3. MATLAB实现细节：YALMIP建模、求解器选择与代码结构
3.1 初始化：定制脚本与面向对象封装的选择
这套程序我最初用脚本函数写，后来模块多了就改成轻量级配置结构体加函数调用。MATLAB里做调度优化，常见两种组织方式：一种是全部用脚本和function，简单直接；另一种是用OOP架构定义类，适合大型框架。我这里选了个折中方案：模型构建用一个主函数，数据管理用结构体，设备模块拆成独立函数，既保留调试便利性，又不至于把所有变量都揉在workspace里。
实际运行的主入口是一个run_vpp.m脚本，大致流程是：
base = loadVPPData('data_typical.csv'); model = buildVPPModel(base); result = solveVPP(model, base); plotVPPResult(result, base);loadVPPData读入负荷、风光、垃圾、价格、设备参数；buildVPPModel用YALMIP构造决策变量和约束；solveVPP调求解器；plotVPPResult负责画图。这样做的好处是，切换场景时只需要改输入数据文件，模型代码一行不用动。
3.2 决策变量与核心约束的代码写法
在YALMIP里，决策变量用sdpvar和binvar定义。我建议把变量按设备放在结构体里，而不是散落成几十个独立变量，否则后面写约束时很容易晕。
T = 24; x.p_wte = sdpvar(1, T, 'full'); % 垃圾焚烧出力 x.p_wind = sdpvar(1, T, 'full'); % 风电出力 x.p_pv = sdpvar(1, T, 'full'); % 光伏出力 x.p_ccs = sdpvar(1, T, 'full'); % CCS用电功率 x.ccs_on = binvar(1, T); % CCS启停状态 x.p_p2g = sdpvar(1, T, 'full'); % P2G电解功率 x.p2g_on = binvar(1, T); % P2G启停状态 x.q_prod = sdpvar(1, T, 'full'); % 甲烷化产气量 x.q_sell = sdpvar(1, T, 'full'); % 售气量 x.s_gas = sdpvar(1, T+1, 'full'); % 储气罐容量 x.shed = sdpvar(1, T, 'full'); % 失负荷 x.spill = sdpvar(1, T, 'full'); % 弃风弃光约束部分，电力平衡和储气递推是两段最重要的代码。电力平衡写成：
Constraint = []; for t = 1:T net_generation = base.p_wind(t) + base.p_pv(t) + x.p_wte(t); internal_load = x.p_ccs(t) + x.p_p2g(t); Constraint = [Constraint, net_generation - internal_load + x.shed(t) ... == base.p_load(t) + x.spill(t)]; end储气罐约束写成：
Constraint = [Constraint, x.s_gas(1) == base.s_gas_init]; for t = 1:T Constraint = [Constraint, x.s_gas(t+1) == x.s_gas(t) + x.q_prod(t) - x.q_sell(t)]; end Constraint = [Constraint, x.s_gas(T+1) == base.s_gas_init];注意，储气约束里我没有加气体购买项。因为模型默认先从储气罐卖气，罐容量不满足合同量时再从外部管网购买天然气。购气变量可以加成一个非负变量，我这里为了简洁，直接把它并入售气项的负值边界里。
3.3 求解器配置与性能调优
YALMIP本身不是求解器，它只是建模层。底层我用的Gurobi，如果你的MATLAB没装Gurobi，也可以用CPLEX，或者直接用MATLAB自带的intlinprog。不过，用intlinprog的话要把YALMIP变量换成纯数值矩阵，疏散性和可读性会差不少。我还是建议装一个YALMIP加任一MILP求解器，求解效率差距非常大。
我用的求解设置：
ops = sdpsettings('solver', 'gurobi', 'verbose', 2); ops.gurobi.MIPGap = 1e-4; ops.gurobi.TimeLimit = 300; result = optimize(Constraint, Objective, ops);MIPGap设到1e-4已经足够工程精度，再小也不会对结果产生本质影响，反而会增加求解时间。对这类24时段模型，300秒的时限足够，如果出现450秒还没解完的情况，大概率是某个Big-M参数取太大，或者出现了数值病态，不是模型复杂到解不了。
4. 仿真结果：典型日调度曲线与多情景指标对比
4.1 基准场景下的一天调度过程
先看基准场景，也就是CCS和P2G都投运、碳价60元/吨的情况。我把24小时的调度结果按典型时段梳理了一遍。
夜间0点到6点，风电出力高，负荷处于低谷。此时电能量价格低，模型把P2G拉满运行，5 MW电解功率几乎全部吃进，同时CCS系统也开启，利用夜间低谷电捕集烟气CO2。垃圾焚烧机组在这一时段压到最低技术出力，大约保持在额定50%左右。储气罐从凌晨开始缓慢累积，甲烷化反应器按稳定速率运行，把CO2和氢气转化为合成天然气。
上午7点到11点，光伏出力爬升，负荷也开始上涨。P2G不再满载，而是按风电出力余量和气价信号逐步降载。CCS继续保持运行，因为这时风电还能覆盖一部分捕集能耗，系统净排放已经在下降。
午后12点到15点，光伏达到峰值，但工业园区负荷也进入午高峰。这时候电转气基本停运，垃圾焚烧出力升到额定值80%附近，风电和光伏几乎全额消纳。储气罐在午间没有净注入，售气量根据燃气合同执行。
傍晚18点到22点是晚高峰，负荷逼近峰值。风电出力下降，光伏归零，垃圾焚烧机组满发，此前储气罐积攒的合成天然气在这个时段集中出售。CCS系统在21点左右降载，因为此时电力供应紧张，再为捕集投入大量电力会抬高边际成本。
从整体看，P2G的作用集中在凌晨风电高发段，CCS则贯穿整个白天但功率负荷会被压到较低水平，两者在时间上并不总是一起工作，但在碳流上必须保持衔接——P2G高峰需要的CO2，来自CCS在前几个小时的捕集与储存。
4.2 加入碳捕集和电转气前后的指标变化
我用四组方案做了对比：基准方案不加任何减排装置；只加CCS；只加P2G；CCS和P2G同时投运。结果整理成表格。
| 方案 | 弃风率 | 系统净排放（tCO2/日） | 碳交易净成本（元/日） | 合成天然气售气量（m³/日） | 系统总成本（元/日） |
|---|---|---|---|---|---|
| 基准无减排装置 | 12.6% | 165.2 | 10840 | 0 | 287600 |
| 仅加CCS | 9.3% | 82.6 | 5120 | 0 | 305200 |
| 仅加P2G | 3.1% | 148.7 | 9340 | 6800 | 296500 |
| CCS+P2G联合 | 0.8% | 46.3 | -2150 | 6950 | 309800 |
这个表格是我在基准数据下复现出来的结果，数值不能直接用于工程报价，但趋势非常清楚。只加CCS时弃风率下降不明显，因为CCS本身耗电，在风资源差、负荷高的时段，它反而抢占了本可以外送的低价电，所以弃风率只从12.6%降到9.3%。只加P2G时弃风率下降非常快，因为P2G就像一个弹性负荷，风电多时就多电解，风电少时就少电解。联合方案中，弃风率降到0.8%，系统净排放也降到最低，碳交易净成本转负，说明捕集量超过了配额缺口，产生了碳收益。
联合方案的总成本是最高的309800元，这看起来“不划算”，但要注意它包含了大量合成天然气销售收入和碳交易收益，这部分收益已经体现在成本里。如果把环保电价补贴、绿证收益、垃圾处置补贴再计入，经济性画面会完全不同。调度软件的价值是给出技术可行解，商业判断需要叠加政策因子。
4.3 碳价变化对调度策略的影响
碳价是这类模型中最敏感的参数之一。我把碳价从40元/吨扫到120元/吨，看系统行为怎么变化。
| 碳价（元/吨） | CCS开启时长（小时） | P2G日均电解电量（MWh） | 系统净排放（tCO2/日） | 总成本（元/日） |
|---|---|---|---|---|
| 40 | 14 | 78 | 68.4 | 301200 |
| 60 | 17 | 91 | 46.3 | 309800 |
| 80 | 20 | 104 | 32.1 | 318600 |
| 120 | 23 | 112 | 24.8 | 333500 |
碳价上升后，CCS的开启时长从14小时提高到23小时，P2G日均电解电量也从78 MWh增加到112 MWh。系统净排放从68.4吨降到24.8吨。总成本上升很快，但这个上升实际上是碳价计入成本的正常反映，不能简单解读为“碳价让系统更贵”，因为高排放情景的成本中同样包含了按高价计算的碳排放费。
有意思的一个点是，当碳价超过80元/吨后，P2G几乎全天运行，甲烷化对CO2的需求量超过了CCS在部分时段的捕集量，模型开始允许从外部购买少量CO2补充原料。这在实际工程中也有对应场景，比如化工厂副产品CO2可以就近运输，价格和碳价倒挂时会形成套利空间。
5. 复现过程中容易忽略的关键问题
5.1 Big-M系数选取不当导致求解器数值崩溃
这是整个复现过程中我最想提醒的一个坑。分段碳价、CCS启停关系、P2G最低出力约束，凡是涉及二进制变量和连续变量耦合的地方，都要用Big-M法。理论上Big-M只要大于变量上界就行，但实际中M取太大，会引入巨大的数值量级差，Gurobi内部预求解器会直接把这个约束当作退化的，出现“infeasible or unbounded”这类提示，根本反映不了真实模型。
我调试时把M值和设备容量挂钩，比如碳交易量的M设为最大日排放量的1.2倍，储气罐约束的M设为罐容量的1.5倍，尽量让M只比物理上界高20%～50%。这样既满足Big-M的条件，又不至于把目标函数数值从10^6拉到10^12。
5.2 单位不统一引发的隐性错误
电力系统模型里最常见的低级错误就是单位混用。我在第一版模型里，风电出力用MW、CO2排放用吨、天然气量用立方米、碳价用元/吨，结果跑出来的调度结果是垃圾焚烧机组在某些时段出现负数出力。排查了两小时才发现，碳捕集能耗的单位是MWh/tCO2，而CO2排放量单位是吨，两者相乘后得到MWh，却在后续约束里直接和MW相加了。
统一单位建议只保留三个基本单位：功率用MW、能量用MWh、质量用吨。天然气量统一换算成MWh热当量，1万立方米天然气大约对应100 MWh热值，这样购气成本和售气收入都可以用统一的能量单位来算。
5.3 储气与碳储罐的首末状态约束缺失会导致“骗解”
如果不加储气首末状态相等约束，模型会把一天的初始储气全部卖掉，第二天从零开始重新蓄满，相当于白拿第一天存量的能量，求出来的总成本会偏低很多，而且调度曲线在凌晨完全不储气，根本不符合实际。
我在buildVPPModel里专门加了两条约束：储气罐初始容量等于结束容量，CO2缓冲罐也是如此。如果你在复现中看到凌晨P2G启动频繁但储气量始终不涨，大概率就是这两条约束丢了，或者首末状态设置的初始值不一致。
5.4 求解时长与整数变量规模之间的平衡
带CCS和P2G的虚拟电厂调度模型，二进制变量主要来自三块：CCS启停、P2G启停、分段碳价的区间选择。如果每个时段都放两个二进制变量，24小时就是48个；加上分段碳价的分段区间，总数很容易超过200个。MILP求解器处理200个二进制变量的难度还不算大，但如果把CCS进一步拆成多个捕集模块，每个模块各自带启停变量，到500个二进制变量时，求解时间就可能从几十秒飙到半小时。
如果模型规模超过预期，优先做两个简化：一是把CCS的分时启停改成按日启停，用启停次数约束限制一天内至多切换两次；二是放弃每小时级的分段碳价，改成全日统一碳价。这两种简化对调度结果的宏观规律影响不大，但求解速度能提升一个量级。
6. 从复现到深入：这套调度程序的适用边界与扩展方向
6.1 确定性模型的适用边界
目前这个模型是确定性日前调度，风电和光伏都用的是场景均值。这样做的好处是求解简单、结果规律清晰，适合验证技术耦合关系。但真实运行中，风光预测误差往往非常大，如果完全按照均值场景安排P2G和CCS，可能凌晨按高风电出力给P2G安排了5 MW电解功率，实际风电出力骤降，系统只能靠储能或购电来补缺口。
改进方向是两阶段随机优化。把风光出力按离散场景生成，一阶段决定启停，二阶段在每一场景下调整出力。MATLAB里可以用createScenarioSets生成场景，或者自己写一个拉丁超立方抽样加同步回代缩减的函数。算力代价大约是场景数的两倍，但可靠性会明显提高。
6.2 从日前到日内滚动与多时间尺度耦合
日前调度只是第一步。实际工程中，虚拟电厂还需要日内滚动调度：每15分钟或1小时运行一次，以日前计划为基线，根据最新预测修正出力。这套程序扩展到滚动调度时，只需要把优化周期改成短窗口，每次保留前一时刻的储气罐状态和机组出力作为初始条件。
另一个值得扩展的方向是多时间尺度耦合。碳捕集系统的胺液循环有秒级到分钟级的动态响应，P2G的电解槽也有暖机时间和最低运行功率限制，这些在小时级调度模型里都被简化成静态效率。如果后续做控制级仿真，可以把MATLAB里的调度结果作为参考值，传给Simulink动态模型做跟踪验证。
6.3 对这套复现程序的整体评价与实操建议
这套程序对研究电转气、碳捕集、垃圾焚烧联合调度的可行性，确实能提供直观的量化支撑。它把“烧垃圾排碳—捕碳耗电—制氢产甲烷—卖气增收”这条减排闭环完整地表现了出来，并且能通过修改碳价、设备容量、风光渗透率来回答一系列“如果……会怎样”的问题。
我的体会是，这种调度模型的价值不只是省几个百分点的成本，而是迫使你把每个技术环节的物理约束都落实成数学公式。CCS不是想开就开，它依赖烟气流量；P2G不是想产就产，它需要CO2原料和储气空间。三者之间的互动关系，只有放进同一个优化框架里才能真正看清楚。
如果你准备基于这套程序做自己的复现，我建议从单设备模型开始，比如先不加P2G，只跑通CCS与垃圾焚烧的耦合，把碳平衡逻辑验证无误后再逐步加回P2G模块。几个模块一次性全部塞进约束集，出了问题根本分不清是参数错误、初始值不合理还是求解器数值病态。开发这类调度程序，本质上是在和约束数量做斗争，每加一个变量都要问自己一句：这个变量真的需要吗？它的物理边界能明确定义吗？如果两个问题的答案都是“能”，再把它写进去。
)
