1. 为什么是"阶梯碳交易+P2G-CCS耦合+燃气掺氢"这个组合
1.1 单一碳价不够用,阶梯碳价才贴近实际碳市场
先说碳交易这部分。看过几篇碳交易相关的论文就会注意到,早期的模型基本都是单一碳价,也就是说,不管你的碳排放量超了多少,每吨碳的惩罚价格是固定不变的。这种处理方式最大的问题在于:它没法体现碳市场"排放越多,边际成本越高"的调控逻辑。现实中的碳配额交易价格是浮动的,碳排放量越大,企业需要购买的额外配额就越多,市场紧张程度上升,价格自然会被推高。所以你看近年EI、SCI期刊上关于低碳调度的论文,很多都在用阶梯碳价(碳价随排放区间递增),目的就是把这种边际惩罚递增的机制写进模型里。
阶梯碳交易的计算逻辑很简单,就是对碳排放量划分几个区间,每个区间对应不同的碳价。比如排放量在配额以内是0成本,超出配额0~500吨的部分按100元/吨,500~1000吨的部分按150元/吨,再往上200元/吨。这样一个分段线性函数,既符合实际碳市场的价格形成机制,又能用混合整数线性规划(MILP)求解,不至于引入非线性让模型难解。
在我复现的这个模型里,阶梯碳交易是整个调度框架的"指挥棒"。它直接决定了虚拟电厂里各设备的出力优先级:碳价越高,碳排放量大的燃气机组就越不划算,系统就会倾向于用风电、光伏这些零碳电源,同时把碳捕集装置开起来。你会发现,碳价的阶梯幅度设置直接影响了调度结果——差值设得太小,碳交易约束形同虚设;设得太大,燃气机组几乎不敢出力,系统又失去了灵活性。这块参数在论文里一般不会给太细的敏感性分析,建议你自己做几组对比实验。
1.2 P2G和CCS耦合:把"碳"变成"甲烷"的闭环逻辑
P2G(Power to Gas)和CCS(Carbon Capture and Storage)单独拎出来都不算新鲜,但把两者耦合起来做虚拟电厂调度,这里面的门道就多了。
先看各自的功能。P2G的核心设备是电解槽和甲烷化反应器:电解槽用多余电能电解水制氢,氢气再和二氧化碳在甲烷化反应器里合成甲烷(也就是天然气的主要成分)。CCS这边是碳捕集装置,从烟气里把CO2分离出来,要么封存,要么送去利用。P2G-CCS耦合的思路,就是让CCS捕集下来的CO2直接作为P2G甲烷化环节的原料,这样既减少了系统的净碳排放,又生产了可存储、可燃烧的天然气,还消纳了原本可能被弃掉的风电。
这个耦合在调度模型里的表现很有意思。传统CCS设备本身是个负荷,因为它要耗电来运行(吸收塔、再生塔、压缩机都要用电),这会降低系统的净出力。但把它和P2G接在一起后,捕集下来的CO2有了去处,甲烷化产生的天然气又能被燃气轮机烧掉,整个能量流和碳流就形成了一个内部闭环。用大白话说,虚拟电厂内部自己组了一个"碳循环工厂"——风电多了就制氢,制出来的氢一部分给燃气轮机掺着烧(这又关联到燃气掺氢),一部分用来和CO2合成甲烷,甲烷存起来等到负荷高峰再发。
在Matlab建模时,这个耦合关系体现在设备之间的能量平衡约束和碳流平衡约束上,CO2捕集量、甲烷化消耗量、封存量三者之间必须满足一个平衡关系式。第一次搭这个约束的时候容易漏掉能量守恒:CCS捕集1吨CO2到底要耗多少电,P2G合成1立方米甲烷要消耗多少CO2和氢气,这些参数都得查文献或者按论文给的数据取值,不能瞎拍。
1.3 燃气掺氢:不是简单把氢混进去烧
燃气掺氢这个技术在近两年的论文里出现频率非常高,背景是天然气管道和燃气机组未来要适应高比例可再生能源系统,氢气作为一种零碳燃料,可以直接按一定比例掺入天然气中使用。对于虚拟电厂场景,掺氢的价值在于给P2G产出的氢气提供了一个就地消纳的出口——氢储能不好搞、储氢罐成本高,但燃气轮机掺氢燃烧就简单得多,把氢气按体积比5%~20%混进天然气里烧,对机组改动不大,却能实打实降低碳排放。
建模的时候,燃气轮机的燃料成本不再是简单的天然气购气成本,而是要同时考虑天然气和氢气两种燃料的消耗量、各自的单价,以及掺氢比例带来的热值变化。另一个关键是燃机出力上限会随掺氢比例变化,因为氢气热值低(按体积算约为天然气的1/3),掺多了会降低机组出力能力。这个细节特别容易被忽略,很多刚上手复现的人直接在约束里写死燃机最大出力,结果算出来的结果跟论文对不上,原因就在这儿。
在这套模型里,掺氢不是固定的,而是可调的决策变量——系统会根据风电出力、氢价、碳价等情况动态决定掺氢比例。风电多的时候多制氢、多掺氢,减少天然气购买;碳价高的时候多掺氢,降低碳排放成本。这就是"优化调度"的灵魂所在。
2. 虚拟电厂整体架构与优化模型构建
2.1 虚拟电厂的设备组成和能量流转关系
这套模型里的虚拟电厂聚合了六类单元:风力发电机组(WT)、光伏机组(PV)、燃气轮机(GT)、电转气设备(P2G,含电解槽和甲烷化反应器)、碳捕集与封存设备(CCS),以及储能装置(包括蓄电池和储氢罐)。它们通过母线连接,对外表现为一个整体,参与电力市场的交易。
能量流方向是这样的:风电和光伏优先给负荷供电,多余电量给蓄电池充电,或者送给P2G电解槽制氢。P2G产出的氢气有两个去向:一是进入储氢罐储存,二是直接送给燃气轮机掺氢燃烧;甲烷化反应器消耗一部分氢气和CCS捕集的CO2合成天然气,送入燃气轮机或者直接卖给气网。燃气轮机发电时产生烟气,烟气进CCS捕集装置,捕下来的CO2回到甲烷化环节,没捕到的部分就是净碳排放,需要购买碳配额。
做到这一步,你就能理解为什么要用虚拟电厂这个框架了——单看P2G或者CCS都是一个独立的设备,但放在虚拟电厂里,它们可以通过协调调度实现比各自独立运行更高的整体经济性和低碳性。负荷低谷期电价便宜,多制氢存起来;负荷高峰期电价贵,燃气轮机掺氢满发,储能放电,多能互补的意义就在这里。
2.2 目标函数:运行成本与碳交易成本的博弈
这个优化模型的目标函数是最小化系统总运行成本,包括购电成本、燃料成本(天然气和氢气)、运维成本、碳交易成本,再扣掉向电网售电的收益。公式用Matlab写的话大概是这样:
objective = sum(购电成本 + 燃料成本 + 运维成本 + 碳交易成本 - 售电收益);这里最值得展开的是碳交易成本的计算方式。系统实际碳排放量减去无偿分配的碳配额,差值就是需要购买的碳配额量。如果实际排放低于配额,多余的部分可以卖出获利。而阶梯碳价的处理方式,是把碳排放区间切成几段,每一段对应一个递增的碳价。
举个例子,假设无偿配额是1000吨,实际排放1600吨,超出600吨。阶梯划分为:超排0~300吨按100元/吨,300~600吨按150元/吨,600吨以上按200元/吨。那么碳交易成本就是300×100 + 300×150 = 75000元。这种分段线性的成本函数直接用线性表达式就能表达,配合0-1变量标记超排量所处的区间,最后交给求解器做MILP求解。
我在写目标函数的时候建议把各项成本分开写成独立的变量或者表达式,不要揉在一起。这么做的好处有两层:一方面方便你调试——结果不对的时候能快速定位到是哪一项成本出了问题;另一方面后面画图分析的时候,你可以很直观地看碳交易成本在总成本中的占比变化,这在论文里是一张很关键的图。
2.3 约束条件:从功率平衡到碳流平衡
约束条件是这个模型最庞杂的部分,我按类别梳理一下:
电力平衡约束:系统内所有电源出力减去设备耗电,加上购电,要恰好等于负荷。注意P2G电解槽和CCS都是耗电设备,燃气轮机的厂用电也要算进去。在Matlab里用等式约束向量化表达:
Constraints = [Constraints, P_wt + P_pv + P_gt + P_dis + P_buy - P_ch - P_elec - P_ccs == P_load];设备出力约束:风电、光伏出力在预测值范围内可调(弃风弃光),燃气轮机出力和爬坡约束,储能充放电约束和SOC递推关系,储氢罐容量约束等。这些约束的注意点是各设备的效率参数要单独定义成参数,方便敏感性分析。
P2G运行约束:电解槽制氢量和耗电量之间通过电解效率耦合;甲烷化环节的耗氢量和产气量、CO2消耗量之间是化学计量比关系。通常假设电解槽制氢效率在70%~80%,甲烷化效率在75%左右,具体数值文献里差别不大。
CCS运行约束:捕集量等于烟气中CO2含量乘以捕集率,捕集过程本身的能耗系数大约是0.2~0.3 MWh/tCO2。CCS的设备容量决定最大捕集量,还有爬坡速率限制。
燃气掺氢约束:燃气轮机实际燃料量是天然气和氢气的混合,需要满足掺氢比例上限约束(一般不超过20%),同时燃机出力和总燃料热值成正比。这块写起来最容易出错,因为涉及体积分数和质量分数的换算,建议统一用能量单位(MWh)来建模,别混用体积、质量、能量三种单位。
碳流平衡约束:CCS捕集的CO2量等于甲烷化消耗量和封存量之和,系统净碳排放量等于实际碳排放量减去捕集量,这个净碳排放量就是后面阶梯碳交易成本的输入变量。
> 注意:约束之间的时序耦合关系很容易错。储能的SOC和储氢罐的储气量都是跨时段的递推变量,务必确保初始时刻的状态给定,否则求解器会给出错误结果却不会报错,排查起来很费劲。3. Matlab实现:从建模到求解的完整流程
3.1 求解框架搭建:Yalmip+Cplex/Gurobi
这套模型本质上是混合整数线性规划问题(MILP),Matlab下最顺手的方案是Yalmip建模工具箱搭配商业求解器Cplex或Gurobi。Yalmip的语法简单,只需要定义变量、写约束、给目标函数,剩下求解交给求解器。很多人在复现这类论文卡壳,问题往往不在模型理解上,而是卡在环境配置——Yalmip版本不兼容、求解器License没配好、或者Matlab版本太新导致Cplex接口崩了。
我个人推荐用Gurobi,因为它的MILP求解速度比Cplex快不少,尤其是这种带有数百个0-1变量的调度模型。安装步骤网上很多,就不重复了,比较关键的几个配置在下面。
% 初始化Yalmip和求解器 addpath(genpath('E:\yalmip')); % 换成你自己的Yalmip路径 addpath(genpath('E:\gurobi')); % 换成你的Gurobi路径 ops = sdpsettings('solver', 'gurobi', 'verbose', 2); ops.gurobi.MIPGap = 0.001; % 设置MIP间隙为0.1% ops.gurobi.TimeLimit = 600; % 最长求解时间600秒MIPGap这个参数很实用。如果模型规模很大,想快速得到一个近似最优解,可以把MIPGap放宽到0.01甚至0.05,求解时间能缩短一个数量级。论文复现阶段建议设小一点,保证解的质量。
3.2 用Yalmip定义优化变量
Yalmip里定义变量用sdpvar和binvar,24小时调度周期,每个设备的出力变量都是一个24×1的向量。关键的变量定义代码如下。
% 时段数 T = 24; % 连续决策变量 P_wt = sdpvar(1, T, 'full'); % 风电出力 P_pv = sdpvar(1, T, 'full'); % 光伏出力 P_gt = sdpvar(1, T, 'full'); % 燃气轮机出力 P_elec = sdpvar(1, T, 'full'); % 电解槽耗电功率 P_ccs = sdpvar(1, T, 'full'); % CCS耗电功率 H_pro = sdpvar(1, T, 'full'); % 电解槽制氢量 H_mix = sdpvar(1, T, 'full'); % 燃气轮机掺氢量 V_p2g = sdpvar(1, T, 'full'); % P2G天然气产量(折算为热值) C_seq = sdpvar(1, T, 'full'); % CO2封存量 C_met = sdpvar(1, T, 'full'); % 甲烷化消耗CO2量 % 储能变量 SOC_bat = sdpvar(1, T+1, 'full'); % 蓄电池SOC(T+1是为了包含初始状态) P_ch = sdpvar(1, T, 'full'); % 充电功率 P_dis = sdpvar(1, T, 'full'); % 放电功率 V_h2 = sdpvar(1, T+1, 'full'); % 储氢罐储氢量 % 0-1变量(用于阶梯碳价分段) u_carbon = binvar(1, 4, 'full'); % 碳排放所处的阶梯区间标记 % 购售电变量 P_buy = sdpvar(1, T, 'full'); % 从电网购电 P_sell = sdpvar(1, T, 'full'); % 向电网售电定义变量的时候有个小习惯要养成:所有功率变量统一用MW,所有能量统一用MWh;时间尺度是1小时,所以功率数值和能量数值在数值上是相等的,这样后续写约束会省很多事。如果系统里混用了kW和MW,即使变量名命名再规范,也很容易在某个约束里写错量纲,导致结果差了好几个数量级。
3.3 阶梯碳交易约束的分段线性化处理
阶梯碳交易是这套模型里最具技巧性的部分,怎么把分段递增的碳价函数写成MILP可解的线性约束,很多第一次接触的人会卡住。思路是用大M法引入0-1变量,让超排量落入对应的排放区间。
假设系统净排放量为E,排放配额为E0,超排量为ΔE = E - E0(≥0)。三个阶梯区间为[0, E1]、[E1, E2]、[E2, +∞),对应碳价λ1 < λ2 < λ3。把ΔE拆成三段:ΔE1、ΔE2、ΔE3,分别表示落在三个区间的排放量。需要满足:
ΔE = ΔE1 + ΔE2 + ΔE3
约束条件是:只有当ΔE1达到上界E1后,ΔE2才能取正值;只有当ΔE2达到上界(E2-E1)后,ΔE3才能取正值。这用线性约束表达就是:
% 分段排放量变量 dE1 = sdpvar(1, 1, 'full'); dE2 = sdpvar(1, 1, 'full'); dE3 = sdpvar(1, 1, 'full'); % u1、u2、u3是0-1变量,标记当前排放处于哪个阶梯 Constraints = [Constraints, dE1 >= 0, dE1 <= E1 * u1]; Constraints = [Constraints, dE2 >= 0, dE2 <= (E2 - E1) * u2]; Constraints = [Constraints, dE3 >= 0, dE3 <= M * u3]; Constraints = [Constraints, dE1 + dE2 + dE3 == E - E0]; Constraints = [Constraints, u1 + u2 + u3 == 1]; % 碳排放成本 C_carbon = lambda1 * dE1 + lambda2 * dE2 + lambda3 * dE3;这个约束的精髓在于:当u2=1时,必须同时有u1=1,也就是前面区间必须达到满载。上面这个写法其实隐藏了这个逻辑——因为U1、U2、U3的和为1,且区间是顺序填满的,通过MIP求解器的分支定界机制能自动找到正确的分段组合。
> 注意:大M的取值很关键。M取太大了会引入数值稳定性问题,取太小了可能限制可行域。在这套模型里,M取800~1000就足够,因为系统最大超排量不会超过这个量级。实际调试时如果发现求解器报数值警告,优先检查大M取值。3.4 燃气掺氢约束的建模细节
燃气轮机的燃料由天然气和氢气混合组成。设天然气热值为9.7 MWh/kNm³,氢气热值为3.0 MWh/kNm³(按体积),掺氢比例为α(体积分数)。燃机输入的混合燃料总热值Q_fuel和出力P_gt之间满足:
Q_fuel = P_gt / η_gt
其中η_gt是燃机效率。混合燃料热值构成关系:
Q_fuel = V_ng × 9.7 + V_h2 × 3.0
掺氢比例约束:
V_h2 / (V_ng + V_h2) ≤ α_max
把这些直接写进Matlab约束:
% 天然气和氢气燃料量 V_ng = sdpvar(1, T, 'full'); V_h2 = sdpvar(1, T, 'full'); Constraints = [Constraints, V_ng * 9.7 + V_h2 * 3.0 == P_gt / eta_gt]; Constraints = [Constraints, V_h2 <= alpha_max * (V_ng + V_h2)]; % 燃料成本 C_fuel = sum(V_ng .* price_ng + V_h2 .* price_h2);这里面有个问题值得留意:氢气来源是P2G,天然气来源是外购或者P2G甲烷化产物。如果甲烷化产物也进入燃气轮机,那V_ng里其实有一部分是系统自产的,成本应该按制取成本算而不是按外购天然气价格算。建模的时候通常会把V_ng拆成V_gas_buy(外购)和V_gas_p2g(P2G产气)两部分,否则成本计算会多算。
3.5 P2G-CCS耦合约束的实现
耦合约束是这套模型区别于普通虚拟电厂模型的核心,分成三个环节来写:
第一个是电解槽环节,氢气产量等于耗电量乘电解效率除以单位制氢电耗。假设电解效率η_elec=0.75,单位氢气低热值为3.0 MWh/kNm³,那么它们的关系是:
% H_pro的单位是kNm³,P_elec的单位是MW H_pro = P_elec * eta_elec / 3.0 * 1000; % 这里按实际单位换算第二个是甲烷化环节,化学反应式是CO2 + 4H2 → CH4 + 2H2O。按体积比,4份氢气生成1份甲烷,消耗1份CO2。所以:
% V_p2g是P2G产甲烷量,C_met是甲烷化消耗CO2量,H_met是甲烷化消耗氢量 C_met = V_p2g; % 化学反应计量比1:1,单位统一为kNm³时 H_met = 4 * V_p2g; % 4:1,需要注意氢气的能量单位换算第三个是CCS环节,CO2捕集量等于燃机烟气中CO2含量乘捕集率,并满足:
% C_cap是捕集量,C_met是甲烷化利用量,C_seq是封存量 C_cap = C_met + C_seq;这三个约束串起来就是P2G-CCS耦合的核心。实际调试时你会发现,CCS捕集率和甲烷化消耗CO2量之间的平衡很敏感,如果捕集量不足或者甲烷化容量受限,多余的风电就只能弃掉或者给蓄电池充电,系统运行经济性会明显下降。这在结果里通常表现为弃风率上升、碳交易成本上升,两个指标一起看就能判断耦合约束是否被正确激活。
4. 案例设置、结果分析与参数敏感性
4.1 算例参数与数据配置
论文复现最烦的一步是找全参数。这类论文一般不会在正文里把所有参数列全,很多藏在附录或者网上公开的数据集里,没有的话就得参考同类文献的取值。我的做法是建一个参数结构体,把所有参数集中在同一个地方定义,这样调整和检查都方便。
%% 参数初始化 para.T = 24; % 调度周期 para.dt = 1; % 时间步长 para.P_load = [/* 负荷数据,1x24 */]; para.P_wt_forecast = [/* 风电预测出力 */]; para.P_pv_forecast = [/* 光伏预测出力 */]; % 设备容量 para.P_gt_max = 200; % 燃气轮机最大出力,MW para.P_gt_min = 20; % 燃气轮机最小技术出力,MW para.ramp_gt = 50; % 燃气轮机爬坡速率,MW/h para.P_elec_max = 100; % 电解槽最大输入功率,MW para.C_ccs_max = 80; % CCS最大捕集量,t/h para.cap_bat = 100; % 蓄电池容量,MWh para.eta_bat = 0.95; % 蓄电池充放电效率 para.cap_h2 = 500; % 储氢罐容量,kNm³ % 价格参数 para.price_buy = [/* 分时购电价,元/MWh */]; para.price_sell = [/* 分时售电价,元/MWh */]; para.price_ng = 2.5; % 天然气价格,元/kNm³ para.price_h2 = 1.8; % 氢气价格,元/kNm³ % 碳交易参数 para.e0 = 500; % 无偿碳配额,t para.lambda = [100, 150, 200]; % 阶梯碳价,元/t para.carbon_interval = [300, 600]; % 阶梯区间边界,t这些参数看起来简单,但每一组数值背后都需要有文献支撑。特别是碳价和配额,不同论文差距很大,有的用欧盟碳市场的价格数据,有的用国内碳试点市场的成交价,数值差个两三倍都很正常。复现的时候务必参考你打算投稿期刊的近期文献,否则审稿人第一个问题可能就是"碳价参数取值依据是什么"。
4.2 求解结果展示与分析切入点
求解完成后,Yalmip会返回目标函数值、变量解和求解器状态。用value()函数提取各变量的数值,然后就可以画图分析了。通常论文里必须有这么几张图:
第一张是电力平衡图,展示24小时内风电、光伏、燃气轮机、储能充放电、购售电的功率分配情况。这张图能直观看到系统的削峰填谷效果,也能验证功率平衡约束有没有写对——把各条曲线求和后应该恰好等于负荷曲线。
第二张是碳流图,展示CCS捕集量、甲烷化CO2消耗量、封存量以及系统净碳排放量的时序变化。这张图直接体现P2G-CCS耦合的运行状态,尤其能看出系统中"碳闭环"是否真正跑通。
第三张是总成本构成饼图或直方图,分解为购电成本、燃料成本、运维成本、碳交易成本、售电收益。通过对比不同类型成本的大小关系,能判断系统低碳运行的主要经济驱动力是来自阶梯碳交易的惩罚还是来自燃料替代的收益。
我复现时的数据结果验证了阶梯碳交易和P2G-CCS耦合的协同效果:碳价阶梯设置越陡,燃气轮机掺氢比例越高,CCS捕集率也越高,但单位减排成本会逐渐上升,存在一个经济最优的碳价水平。这就是论文里经常讨论的"碳交易机制设计对调度的影响"的切入点。
4.3 敏感性分析怎么做才够论文标准
审稿人对纯单算例的结果通常不满意,需要补充敏感性分析来证明模型规律的普遍性。这套系统里最值得做敏感性分析的就三个参数:碳价、掺氢比例上限、P2G设备容量。
碳价敏感性:把阶梯碳价从50/80/120逐档提升到300/400/500,观察系统碳排放总量和总成本的变化趋势。一般在图纸上会看到碳排放总量随碳价上升而单调下降,但下降速率越来越慢,呈现边际递减效应,这就说明存在一个合理的碳价区间,超过这个区间继续提高碳价对减排的促进效果有限,反而大幅增加运行成本。
掺氢比例上限敏感性:从5%到20%以5%为步长扫描,观察燃气轮机出力变化、P2G设备利用率和系统碳排放。掺氢比例上限提高后,系统能消纳更多氢气,P2G的利用率提升,弃风率下降。但如果储氢容量有限或者掺氢比例过高导致燃机效率下降,总成本可能出现拐点。
P2G容量敏感性:从50MW扫到200MW,观察系统总成本和碳排放的变化。P2G容量增大意味着能消纳更多廉价风电,但同时设备投资成本上升,折合到日运行成本里的折旧费用也更高等。
> 实操心得:敏感性分析一定要写成循环批量跑,不要手动改参数一次次运行。外层循环写参数值,内层调用同一个求解函数,结果存到结构体里统一画图。一个参数5个档位,每组求解一两分钟,总共十几分钟就跑完了。5. 常见问题与排查技巧实录
5.1 模型无解或者求解器报"infeasible"
这是复现阶段最常碰到的问题。模型无解的原因按出现频率排,大致是这几种:
第一类是约束自相矛盾。最典型的是把储能SOC初始状态和末尾状态同时设为固定值,但充放电功率上限又不足以在调度周期内完成充放电切换,导致可行域为空。再比如燃气轮机的爬坡约束设得太紧,而负荷曲线波动大,系统无法跟上负荷变化。排查方法很简单:加约束前先注释掉一部分,逐组排查哪个约束组导致无解。在Yalmip里用check(Constraints)命令查看各约束的残差,能精确定位到具体的约束行。
第二类是数据问题。负荷或新能源预测数据里有异常值(比如某时刻负荷为负、风电出力超过装机容量),这些数据会直接让模型无解。写代码前先做数据清洗,plot一遍所有输入数据看有没有明显离谱的点。
第三类是参数量纲混用,这个前面强调过。kW和MW混用、MWh和MJ混用在数值上差了10⁶,约束基本不可能有解。用求解器之前花两分钟把参数的量纲统一,能省一晚上调试时间。
5.2 求解速度慢,MILP分支定界树爆炸
系统规模不大(单时段24个点、设备几十个),理论上MILP求解应该很快。如果求解速度明显慢,通常是这几个原因:
0-1变量太多。阶梯碳交易的分段如果设得太多,比如5段以上,每个状态都要引入0-1变量,分支定界速度会明显下降。建议碳排放区间分段控制在3~4段,既能体现阶梯机制,又不会让求解太慢。
约束冗余。比如给每个时段都写了全时段的储能SOC约束,但有些约束在数学上是重复的,会导致求解器预处理时做太多冗余运算。检查有没有重复约束,把重复的Constraints = [Constraints, ...]删掉。
求解器参数没调。MIPGap从默认的1e-4放宽到1e-3,求解时间可能差十倍。TimeLimit设一个上限,求解到时间上限就取当前最好的可行解,这在论文复现中完全够用。
5.3 结果和论文对不上,差在哪里
这是复现论文的终极难题。我花了很长时间才发现几个常见差异来源:
第一种是碳配额分配方式不同。有的论文用历史法(基于历史排放量按比例分配),有的用基准线法(基于发电量乘以基准排放因子),两种方法算出的配额量差很多,直接导致碳交易成本变化。如果复现的数值对不上,先检查论文是用的哪种配额分配方法。
第二种是碳排放核算边界不同。有的模型只计算燃气轮机燃烧产生的碳排放,有的还计入外购电力的间接碳排放,还有的从整个系统生命周期角度算碳排放。边界不同,碳成本差异非常大。先把论文的碳核算边界看清楚,再去调参数。
第三种是电价机制不同。有的算例用分时电价,有的用实时电价,还有的用固定电价+平衡市场机制。电价曲线不同,储能和P2G的运行策略会完全不同。如果论文给的是某实际电力市场的电价数据,尽量找到原始数据。
> 独家技巧:不要只看论文给的最终数值去对比,要看趋势、看模式。比如"储能大概率在凌晨充电、晚高峰放电""P2G在风电大发时段制氢""碳价提高后燃机出力下降"——先确认这些定性的运行规律复现出来了,再纠结具体的数字差异。运行规律对了,参数微调就能把数字也对上。5.4 画图与数据分析的实用建议
最后说点画图的事。Matlab出图推荐用gcf、subplot和yyaxis的组合,一套图能看出所有信息。我自己习惯的排版是:一行电力平衡图,一行碳平衡图,一行成本构成图,每张图用不同颜色区分不同设备,图例放在外侧避免遮挡数据。
一个提升论文质量的小细节:除了画时序曲线,建议把"有/无阶梯碳交易""有/无P2G-CCS耦合""有/无燃气掺氢"这几种模式下的运行结果做对比,用一张表列出总成本、碳排放、弃风率、掺氢比例这几个关键指标。这就是消融实验的思想,也是这类论文里最能体现贡献的部分,审稿人很吃这一套。
我个人在实际操作中的体会是:复现这类模型,最难的不是读懂数学表达式,而是把论文里说的每一句话精确翻译成代码,哪怕是一个"效率系数取0.85"这样的小参数,都可能让结果相差甚远。建议每篇论文先搭建一个最小可行模型——只保留燃气轮机和风电,跑通之后再逐步加上储能、P2G、CCS、碳交易,每加一个模块就调试验证一次。这种方法虽然前期慢,但后面整体联调时反而省时间,因为问题已经被局部化到最近加入的模块了。说穿了,复现文章拼的不是聪明,而是耐心和系统化debug的工程习惯。
这套模型后续还可以往几个方向扩展:比如加入需求响应机制,让负荷侧也参与调度;或者把日前市场和实时市场两级优化结合进来;再或者考虑不确定性,用鲁棒优化替代确定性调度。不管往哪个方向走,这个"阶梯碳交易+P2G-CCS+掺氢"的基础结构都够用,改起来也顺手,这也是我推荐先花时间把它吃透的原因。