阶梯型碳交易与电制氢综合能源系统MILP优化调度

简介：面向能源系统优化方向的研究人员与具备MATLAB基础的工程技术人员，这份PDF资料围绕阶梯型碳交易与电制氢综合能源系统优化调度展开，重点处理传统机组、风电、电解槽、储氢罐和燃料电池五类单元的协同建模。内容涵盖阶梯碳成本分段线性化、氢能系统状态方程、目标函数设计及YALMIP调用求解思路，并通过算例展示风电大发时段优先制氢、碳排放接近阈值时火电出力调整等现象，附有代码实现与可视化结果，便于读者理解混合整数线性规划在电-氢-电多向转换中的落地方式。资源包仅含1个PDF文件，大小约129KB，结构紧凑，适合快速查阅模型框架与关键约束。已有163人学习，可为综合能源系统调度建模、碳交易机制嵌入及氢储能动态约束处理提供参考。

1. 需求先立住：阶梯型碳交易叠加电制氢，调度问题到底难在哪

做综合能源系统优化调度的人，早晚会撞上一个尴尬场面：模型里碳排放按固定单价算，求解结果在纸面上很漂亮，可政策一旦换成阶梯型碳交易——排得越多单价越贵——原来那条调度曲线立刻失效。再叠加电制氢，电解槽把富余风电转成氢存进储氢罐，燃料电池再把它变回电，电-氢-电这条链路一进模型，变量维度、时序耦合、效率折算方向全变了，MILP 的规模也随之上去。

这套用 MATLAB + YALMIP + CPLEX 写的调度程序，就是把这几个机制塞进同一个混合整数线性规划里，跑 24 小时逐时出力。它适合两类人：正在做能源方向毕设、课设，需要一套结构完整、能被追问的模型；以及要给园区或微网做一版“算得动、说得清”的调度方案的技术人员。参数表、约束构建、阶梯碳价分段、储氢状态方程、求解器调参、结果校验，一条线走完。

2. 系统边界与决策变量：五个单元怎么落到 24 小时矩阵上

模型边界先划清楚，后面所有纠结都源于这一步含糊。这套系统的物理构成是传统火电机组、风电、电解槽、储氢罐、燃料电池五个单元，加上一个外部电负荷和一路氢气负荷。优化周期取 24 个时段，步长 1 小时，不做机组组合的多日滚动，只做日前调度。这个粒度是国内毕设和多数园区方案的通行选择——再细到 15 分钟，风电预测误差会盖过优化带来的收益；再粗到 4 小时，储氢罐的时序耦合基本失效。

2.1 决策变量表与维度约定

先把变量摆成一张表，别急着写代码。变量一旦定下来，约束就是照着表填空。

维度写法上有个习惯：时间维放第一维，设备维放第二维。这样sum(P_gt, 1)拿到的是每台机组的日发电量，sum(P_gt, 2)拿到的是每个时段的系统总出力，写约束的时候不容易搞反。我见过有同学把维度定反，结果爬坡约束变成了“跨机组爬坡”，模型能跑通但物理上完全说不通，答辩时被一句“你这台机组怎么在两台机组之间爬坡”问住。

2.2 参数初始化与量纲统一

参数文件单独放一个sys_init.m，别散在主脚本里。散着写的后果是改一个碳价要翻三个文件，调到第五版自己都不记得哪个数是当前的。

% sys_init.m —— 系统参数集中管理 % 量纲约定：功率 MW，时间 h，氢 kg，碳 t，成本 元 units.dt = 1; % 调度步长 h units.T = 24; % 调度时段数 units.gen_cost = [28 31 25]; % 三台火电机组燃料成本系数 元/MWh units.gen_pmax = [80 60 50]; % 机组出力上限 MW units.gen_pmin = [20 15 12]; % 机组出力下限 MW units.gen_ramp = [30 25 20]; % 机组爬坡限值 MW/h units.wind_cap = 150; % 风电装机容量 MW units.load = 180*ones(1,24) + 40*sin((1:24)/24*2*pi); % 电负荷 MW units.wind_avail = 120*ones(1,24); % 风电可发功率 MW（换成实测风速文件） units.ele_pmax = 20; % 电解槽额定功率 MW units.ele_pmin_r= 0.3; % 最小负载率，低于此值效率骤降 units.ele_eff = 0.72; % 电解槽额定工况效率 units.h2_per_mwh= 20; % 单位换算：理论产氢 20 kg/MWh units.fc_pmax = 10; % 燃料电池发电上限 MW units.fuel_cell_eff = 0.58; % 燃料电池发电效率 units.h2_cap = 2000; % 储氢罐容量 kg units.h2_init = 400; % 初始储氢量 kg units.h2_demand = 120*ones(1,24); % 氢气负荷 kg/h units.ef_gt = 0.85; % 火电碳排放因子 t/MWh units.carbon_seg = [200 500 800]; % 阶梯碳价分界点 t units.carbon_price = [80 120 160]; % 对应档位碳价 元/t units.wind_penalty = 150; % 弃风惩罚 元/MWh

这张参数表里最容易翻车的是h2_per_mwh和carbon_price的量纲。碳价在政策文件里可能写成“元/吨”，也可能写成“元/kg”，两者差 1000 倍；氢气流量在文献里常给“Nm³/h”，转成 kg 还要除以 11.2。我一般会在参数文件顶部写死一行量纲注释，之后每加一个参数就对着这行核一遍，比事后查错误省事得多。

另外wind_cap和wind_avail是两回事：前者是装机容量，用来卡上界；后者是逐时可发功率，来自风速曲线折算。有的代码把两者混成一个变量，结果风电出力曲线是一条直线，优化出来的“消纳策略”全是假的。

2.3 功率平衡与时序耦合的边界条件

功率平衡是整个模型的骨架，写错了后面全废。系统级的电功率平衡是：

% build_model.m —— 功率平衡与风电子平衡 cons = []; % 1) 系统电功率平衡：源 = 荷 for t = 1:units.T cons = [cons, sum(P_gt(t,:)) + P_wind_use(t) + P_fc(t) ... == units.load(t) + P_ele(t)]; end % 2) 风电分解：可发功率 = 上网 + 弃风 for t = 1:units.T cons = [cons, P_wind_use(t) + P_curt(t) == units.wind_avail(t)]; cons = [cons, 0 <= P_wind_use(t) <= units.wind_cap]; cons = [cons, P_curt(t) >= 0]; end

第一组约束里，P_ele(t)挂在负荷侧而不是源侧，因为从电网络看电解槽就是个可调负荷。这一点在写代码时容易想反，写成+ P_ele会让模型认为电解槽在发电，求解器倒是能给出“最优解”，只是那个解物理上不成立——储能越放越多，碳排越低，看着特别美好。

火电机组的上下限和爬坡要分开写，别合成一个约束：

% 机组出力上下限 + 爬坡约束 for t = 1:units.T cons = [cons, units.gen_pmin <= P_gt(t,:) <= units.gen_pmax]; end for t = 2:units.T cons = [cons, -units.gen_ramp <= P_gt(t,:) - P_gt(t-1,:) ... <= units.gen_ramp]; end

首末时段的边界条件必须显式给。储氢罐的日周期闭合是最典型的一处：如果不要求H2(24)回到初值附近，优化器会心安理得地把罐子抽空，用“透支明天的氢”换今天的低成本，日成本确实降了，但这套调度没法连续执行。我一般写成允许小幅偏差的软闭合，偏差量进目标函数罚一笔：

% 储氢日周期闭合（软约束，允许 ±5% 偏差） cons = [cons, H2(24) >= 0.95*units.h2_init, ... H2(24) <= 1.05*units.h2_init];

机组组合类的边界同理，P_gt(1,:)的初始状态要跟上一日末状态对齐，否则第一天永远从“冷态”起算，日成本曲线会呈现一个假的锯齿。

3. 阶梯碳交易的分段线性化：为什么价格递增时连二进制都不用加

阶梯碳价是这个模型里最容易被写歪的一块。它的本质是一条分段线性的凸函数：排放量落在第一档按低价，超出部分按高价，边际成本随排放量跳涨。这条曲线在数学上很简单，在 YALMIP 里却极容易写成非线性。

3.1 直接写 min/max 会毁掉 MILP 结构

很多教程里给的写法是这样的：

% 错误示范：数学上对，YALMIP 里不行 carbon_cost = 80*min(E_total, 200) + 120*max(0, E_total - 200);

如果E_total是一个数值，这行没问题；但它是一个 sdpvar，min和max作用在决策变量上会生成非光滑的分段表达式，YALMIP 要么直接报“nonconvex”，要么把它转成需要非线性求解器的形式。一旦落到 MINLP，CPLEX 就不能用了，得换 BARON 或 SCIP，求解时间从几秒涨到几十分钟，24 时段、三台机组的规模都可能解不动。

更隐蔽的一种写法是用if-else做碳价判断：

% 错误示范二：条件分支对 sdpvar 无效 if E_total <= 200 carbon_cost = 80*E_total; else carbon_cost = 80*200 + 120*(E_total-200); end

这段代码不会报错，但E_total <= 200这个比较对 sdpvar 永远返回 false（或者被当成未知），程序会径直走进 else 分支，把整个模型锁死在第二档。求解能跑完，结果看着也像那么回事，只有在排放量本该落在第一档的算例里才会暴露——曲线对不上，成本偏高，查半天查不到原因。这个坑我踩过，血泪经验是：凡是if后面跟着 sdpvar，一律当成 bug。

3.2 三段辅助变量的标准写法

正确做法是引入分段辅助变量，把非光滑函数摊平成线性约束。设排放量分三档，分界点 200 / 500 / 800 吨，对应碳价 80 / 120 / 160 元每吨：

% 阶梯碳交易：分段线性化 K = numel(units.carbon_price); % 档位数 3 delta = sdpvar(units.T, K, 'full'); % 各时段各档排放量 t % 每个时段的排放量按档位拆开 for t = 1:units.T cons = [cons, sum(delta(t,:)) == units.ef_gt * sum(P_gt(t,:))]; % 第 1 档容量：分界点本身 cons = [cons, 0 <= delta(t,1) <= units.carbon_seg(1)]; % 第 2、3 档容量：相邻分界点之差 for k = 2:K-1 cons = [cons, 0 <= delta(t,k) <= units.carbon_seg(k)-units.carbon_seg(k-1)]; end cons = [cons, delta(t,K) >= 0]; % 最高档不设上限 end % 碳成本 = 各档排放量 × 对应碳价 carbon_cost = 0; for k = 1:K carbon_cost = carbon_cost + units.carbon_price(k) * sum(delta(:,k)); end

关键点在delta的排列顺序上：delta(t,1)装第一档的排放量，容量 200 吨；delta(t,2)装第二档，容量 300 吨；delta(t,3)装第三档，不封顶。约束只有三条——总量等于实际排放、每档非负、每档不超过该档容量。没有任何二进制变量，也没有 big-M。

这样写能成立的逻辑是：碳价 80 < 120 < 160 单调递增，目标函数是最小化总成本，求解器一定优先把便宜的档位填满，再动贵的档位。换句话说，全序关系天然由目标函数保证，不需要额外的整数变量去“选择”档位。这一点很多人没想透，见到分段就上 big-M，MILP 规模凭空翻三倍，求解时间从 3 秒变 3 分钟。

提示：这个“免二进制”的结论只在段价格单调递增时成立。如果某一段是固定费用（比如超过阈值先收一笔基础费），或者存在阶梯补贴导致价格递减，LP 松弛会失效，必须补二进制变量，下面这一节就是那种情况。

3.3 价格非单调时才需要的二进制写法

有些园区的碳配额政策里，超出配额的部分先给一段低价奖励，再进入高价惩罚，价格序列不是单调的。这种模型必须用二进制变量锁定区间，配套 big-M 处理。

% 非单调阶梯：二进制选档 + 动态 big-M z = binvar(units.T, K, 'full'); % z(t,k)=1 表示时段 t 落在第 k 档 for t = 1:units.T cons = [cons, sum(z(t,:)) == 1]; % 每时段只能选一档 E_t = units.ef_gt * sum(P_gt(t,:)); for k = 1:K % 上界直接用该档的物理最大值，而不是拍脑袋的 1e6 M_k = units.carbon_seg(min(k, numel(units.carbon_seg))); cons = [cons, E_t <= M_k + (1 - z(t,k)) * M_k]; cons = [cons, E_t >= M_k * z(t,k)]; end end

这里 M 的取法是整个模型能否稳定收敛的关键。手写1e6是最常见的做法，也是最常见的翻车源头：CPLEX 的分支定界对约束矩阵的条件数敏感，M 比实际量级大三个数量级以上，LP 松弛会变得极松，节点数暴涨，跑一小时还在 gap 20% 上晃。我的做法是用该档的物理上界M_k做动态 M，量级和决策变量对齐，通常把 gap 收到 1% 以内的时间能缩短一个数量级。

4. 电制氢链路：效率双向折算和状态方程里的索引对齐

氢能这块的复杂度不在变量多少，而在效率折算的方向。电解槽是“电进氢出”，燃料电池是“氢进电出”，两个方向在同一个储能方程里相遇，符号或者除数写反一个，模型就会出现负储量这种反物理的解。

4.1 电解槽效率的分段线性化

电解槽的效率不是常数。低负载率下电化学反应效率下降明显，额定负载附近最优，过载又略微回落。原始代码里给了三个工作点：负载率 0.3 时效率 0.6，0.5 时 0.72，1.0 时 0.68。用这三组点做分段线性，比强行用一个常数效率更接近实测曲线。

麻烦在于产氢量是效率乘功率，直接写eff(P)*P是双线性项。标准处理是引入权重变量 λ，对功率点和产氢点同时做凸组合：

% 电解槽分段线性：SOS2 权重法 lam = sdpvar(units.T, 3, 'full'); % 三个工作点的权重 P_pts = units.ele_pmin_r * units.ele_pmax * [1, 5/3, 10/3]; % 6 / 10 / 20 MW eff_pts = [0.60, 0.72, 0.68]; % 各工作点对应的产氢率 kg/h（效率 × 理论产氢 × 功率） h_pts = eff_pts .* units.h2_per_mwh .* P_pts; % 72 / 144 / 272 kg/h for t = 1:units.T cons = [cons, sum(lam(t,:)) == u_ele(t)]; % 停机时权重全零 cons = [cons, P_ele(t) == sum(lam(t,:) .* P_pts)]; % 功率凸组合 cons = [cons, H2_prod(t) == sum(lam(t,:) .* h_pts)]; % 产氢量凸组合 end cons = [cons, sos2(lam')]; % 相邻两点之间插值

sum(lam(t,:)) == u_ele(t)这一行把启停和出力绑在一起：u_ele(t)=0时三个权重必须全零，电解槽功率自然为零；u_ele(t)=1时权重和为一，功率是三个工作点的凸组合。sos2约束保证只有相邻两个权重非零，等价于“落在某一段内插值”，避免求解器把三个点混着用——不加 SOS2，模型会取三个点的加权平均，得到一条并不存在的效率曲线。

参数上面的P_pts是从最小负载率折算出来的：0.3×20=6 MW，0.5×20=10 MW，1.0×20=20 MW。第三点写 10/3 是 20/6 的简化，实际写代码时直接写 20 更清楚，我这里保留系数是为了让你看清折算关系。h_pts的三个值分别是 0.6×20×6=72、0.72×20×10=144、0.68×20×20=272 kg/h，量级和后面的储氢罐容量要对得上。

4.2 储氢罐状态方程与时间索引对齐

储氢罐是标准的储能模型，难点全在索引上：

% 储氢罐动态：时段末储量 = 时段初储量 + 制氢 - 耗氢 - 氢负荷 for t = 2:units.T cons = [cons, H2(t) == H2(t-1) ... + H2_prod(t-1) ... - P_fc(t-1) * units.h2_per_mwh / units.fuel_cell_eff ... - units.h2_demand(t-1)]; end % 容量上下限，留 10% 安全余量 cons = [cons, 0.1*units.h2_cap <= H2 <= 0.9*units.h2_cap];

注意燃料电池这一项的除号位置：P_fc是发电功率，要换算回耗氢量，得先乘单位换算系数再除以发电效率，即P_fc × 20 / 0.58。10 MW 出力对应约 345 kg/h 的耗氢量。如果写成乘效率，耗氢量会变成 116 kg/h，少算了三分之二，模型会以为氢气取之不尽，最优解里燃料电池几乎全天满发——碳成本是降了，但那个解没法执行。

索引对齐是第二个隐雷。储能方程里H2(t)用的是H2_prod(t-1)和P_fc(t-1)，而不是同一时段的功率。这代表“上一时段充进来的氢，在下一时段初到账”。方差一小时，整个日周期的氢平衡会漂移；循环只写到t=24，最后一个时段的制氢量就永远进不了罐子，而它照样烧了电、算了成本。代码注释里那句“注意储氢状态从 t=2 开始，但电解槽功率从 t=1 开始”，说的就是这件事。稳妥的写法是在循环结束后补一条闭合校验：

% 全周期氢平衡校验：制氢总量 = 耗氢 + 氢负荷 + 储量变化 cons = [cons, sum(H2_prod) == ... sum(P_fc) * units.h2_per_mwh / units.fuel_cell_eff + ... sum(units.h2_demand) + (H2(units.T) - units.h2_init)];

这条约束在数学上冗余，但它是一个很好的自检：如果模型能解而这条约束不满足，说明前面的时序索引一定错位了。

4.3 启停次数约束与燃料电池的爬坡

电解槽的启停损耗是真实成本，不能忽略。24 小时内最多启停 3 次，这个约束用状态变量的一阶差分来写：

% 电解槽启停次数约束（24 小时内不超过 3 次） startup = sdpvar(units.T, 1, 'full'); for t = 2:units.T cons = [cons, startup(t) >= u_ele(t) - u_ele(t-1)]; cons = [cons, startup(t) >= 0]; end cons = [cons, sum(startup) <= 3]; % 最小负载率耦合：运行时功率不得低于 30% 额定 cons = [cons, P_ele >= units.ele_pmin_r * units.ele_pmax * u_ele]; cons = [cons, P_ele <= units.ele_pmax * u_ele];

startup(t) >= u_ele(t) - u_ele(t-1)是一个方向约束，只在“从停到开”时右侧为 1，其余情况为 0 或负。因为目标函数里通常会给启停加一个小的惩罚系数（哪怕只有 5 元/次），求解器就不会让startup随意取大值。不加这个惩罚项，startup会退化成无界变量，启停次数约束形同虚设。

燃料电池这边要补爬坡和最小出力，逻辑和火电一样，只是上限更小：

cons = [cons, 0 <= P_fc <= units.fc_pmax * u_fc]; for t = 2:units.T cons = [cons, -8 <= P_fc(t) - P_fc(t-1) <= 8]; % 爬坡 8 MW/h end

最后把所有成本拼进目标函数。这里有个顺序问题：碳成本必须放在循环外面，用全周期的delta求和，不能每个时段单独算一次再累加——那样每段都会从第一档起算，等于把一次阶梯计价拆成了 24 次，碳价被严重低估。

% 目标函数：燃料 + 弃风惩罚 + 阶梯碳 + 氢能运维 + 启停 fuel_cost = sum(sum(P_gt .* units.gen_cost)) * units.dt; curt_cost = units.wind_penalty * sum(P_curt) * units.dt; h2_op_cost = 3 * sum(P_ele) * units.dt + 5 * sum(P_fc) * units.dt; start_cost = 20 * sum(startup) + 20 * sum(startup_fc); total_cost = fuel_cost + curt_cost + carbon_cost + h2_op_cost + start_cost; ops = sdpsettings('solver','cplex','verbose',0, ... 'cplex.timelimit', 3600, ... 'cplex.mip.tolerances.mipgap', 0.01); diagn = optimize(cons, total_cost, ops);

mipgap设 1% 是工程上的通行选择：24 时段调度的时间分辨率就是 1 小时，功率精度到小数点后三位对实际执行没有任何意义，而 gap 从 1% 收到 0.1% 往往要多花几倍时间。timelimit给 3600 秒是兜底，正常算例十几秒就收敛了，设这个值是为了防止某个参数组合把求解器拖死时无人值守。

5. 避坑与排查：五类让模型“能跑但不对”的坑

比起报错，更麻烦的是模型跑通、结果有图、结论却不对。下面这五条是这类项目里出现频率最高的，每条按现象、原因、解决三段说。

5.1 求解状态类：infeasible 和 unbounded 的排查顺序

现象：optimize返回diagn.problem = 3（infeasible）或 4（unbounded），但约束清单看上去没有矛盾。

原因：绝大多数是功率平衡和储能闭合这两组约束打架。最常见的一种是日周期闭合写成硬约束H2(24) == h2_init，而氢负荷已经把罐子抽到下限以下，无论怎么调度都回不到初值。第二种是火电最小出力之和大于负荷峰值，系统“关不掉”机组，多发的电无处可去——因为模型里没有弃电变量。

解决：YALMIP 有个很好用的诊断函数，出问题第一时间跑它，而不是靠肉眼翻约束。

% 不可行诊断：定位到具体哪一组约束冲突 if diagn.problem == 3 [c, idx] = check(cons); disp(c); disp(idx); % idx 给出冲突约束在 cons 里的序号 end

先看idx指向哪几条，再判断是“物理上本来就不可能”还是“参数写错了”。前者要改模型结构（加弃电变量、把闭合约束软化），后者改参数。另一个习惯是把模型规模先降到 4 个时段跑一遍，小规模可行、大规模不可行，八成是爬坡或者容量约束在跨时段上出了问题。

5.2 数值类：big-M 与量纲带来的隐性失真

现象：模型能解，但 MIP gap 长时间停在 15% 以上；或者换个碳价参数，最优解突然大幅跳变，和物理直觉不符。

原因：两种可能。一是 big-M 用了1e6这种量级，约束矩阵条件数恶化，LP 松弛极松；二是量纲不统一，比如功率用 kW、碳价用元/吨、时间用分钟，混在一起之后各成本项相差几个数量级，求解器的尺度变换失效。

解决：big-M 一律用物理上界动态计算，别手写常数。量纲在参数文件顶部统一注释，全部落到 MW / h / kg / t / 元这一套。另外打开verbose=2观察求解日志的前几行，YALMIP 会打印预处理后的变量数和约束数，如果约束数比你预想的多出好几倍，说明有约束被重复生成了——循环里cons = [cons, ...]写重复是常见原因。

% 求解前先做一次规模体检 disp(['变量数: ', num2str(length(getvariables(cons)))]); disp(['约束数: ', num2str(length(cons))]);

5.3 物理类：负储氢和效率倒挂

现象：结果图里储氢量出现负值，或者电解槽在风电零出力的时段满负荷运行。

原因：负储氢说明状态方程的方向写反了，或者容量下界漏了H2 >= 0这一条（只写下限0.1*cap是不够的，某些写法里会被覆盖）。电解槽在无风时段猛开，通常是功率平衡的符号写反，把P_ele加到了源侧，系统凭空多出一个“发电设备”。

解决：加一条“物理校验”脚本，在optimize之后逐条核对基本守恒：

% plot_results.m 里的自检段 H2_val = value(H2); Pe_val = value(P_ele); Pc_val = value(P_fc); assert(all(H2_val >= -1e-6), '出现负储氢，检查状态方程方向'); assert(all(Pe_val <= units.ele_pmax + 1e-6), '电解槽超额定，检查上限约束'); bal = sum(value(P_gt),2) + value(P_wind_use) + Pc_val ... - units.load' - Pe_val; assert(max(abs(bal)) < 1e-4, '功率不平衡，检查平衡式符号');

assert在正式代码里通常被去掉，但调试阶段我强烈建议留着。求解器不会告诉你“这个解不物理”，它只保证约束满足，约束本身写错它一样给“最优解”。

5.4 碳成本记账类：重复起算与循环内累加

现象：碳成本占总成本的比例明显偏低，或者把碳价从 80 调到 200，总成本几乎不变。

原因：碳成本写在了时段循环内，每个时段都从第一档起算，等于给了 24 次低价额度；或者E_total只取了最后一个时段的排放量，前面 23 个时段的排放凭空消失。

解决：碳成本一定用全周期的delta求和，sum(delta(:,k))而不是delta(t,k)。核对办法很简单，把最优解里的碳成本手动重算一遍，和value(carbon_cost)比对：

% 手工复核碳成本 E_manual = units.ef_gt * sum(sum(value(P_gt))); c_manual = 0; for k = 1:numel(units.carbon_price) c_manual = c_manual + units.carbon_price(k) * sum(value(delta(:,k))); end fprintf('模型碳成本 %.2f，手工复核 %.2f\n', value(carbon_cost), c_manual);

两个数对不上，就说明delta的定义或者求和维度有问题。

5.5 结果解读类：把算例结论当成普遍规律

现象：算例跑出来“阶梯碳价下总成本比固定碳价低 8.6%”，写进报告时当成一条普适结论。

原因：这个结论高度依赖参数。原始算例里固定碳价取的是一档偏高的统一价，阶梯机制下系统把大部分排放压在第一档低价区间，所以总成本反而低。如果统一碳价取第一档的 80 元/吨，结论立刻反过来。

解决：汇报时把“参数敏感性”作为结论的一部分。至少扫三组碳价，画出总成本和碳排量的权衡曲线：

% 碳价敏感性扫描 price_list = [60 80 100 120 150 200]; for i = 1:numel(price_list) units.carbon_price = [price_list(i) 1.5*price_list(i) 2*price_list(i)]; build_model; % 重建模型 optimize(cons, total_cost, ops); res.cost(i) = value(total_cost); res.emis(i) = units.ef_gt * sum(sum(value(P_gt))); end plot(res.emis, res.cost, 'o-'); xlabel('碳排放 t'); ylabel('总成本 元');

有了这条曲线，答辩时被问“换个碳价你的结论还成立吗”，直接甩图比解释十分钟管用。原始算例里那两个现象——凌晨风电高峰优先制氢而非弃风、碳排接近阈值时火电出力出现拐点——都能从这条曲线上找到对应，说明模型的行为是可解释的，不是数值噪声。

6. 验收与答辩：怎么证明这套调度结果可信

模型跑出图不算完，能证明这张图可信才算完。验收环节我一般分三层：约束层看可行性，物理层看是不是反常识，经济层看成本项的相对关系。

约束层的验收最机械，也最容易被跳过。把功率平衡残差、储氢上下限、启停次数、碳排分段容量这几条逐项打印出来，误差阈值写 1e-4，全部通过才算模型收敛到了可行域。物理层要盯三张图：各单元 24 小时出力堆叠图，看源荷是否匹配、有没有出现无风时段电解槽满发的怪象；储氢罐储量曲线，看日周期是否闭合、有没有贴死上下限（贴死说明容量或负荷配比不合理）；碳排逐时曲线叠加阶梯分界线，看阈值附近的出力拐点是否真实存在。经济层最简单，看成本构成饼图——碳成本占比通常在 15% 到 20% 之间，低于 10% 说明碳价或者排放因子量纲错了，高于 30% 说明火电配比过高，氢能没起到替代作用。

答辩时被追问最多的两个问题，我提前给个应对思路。第一个是“为什么不直接上非线性求解器，把电解槽效率建成二次函数”，回答要落到可解性上：24 时段、三机组、含二进制的 MINLP，CPLEX 用不了，换非线性求解器不仅求解时间不可控，还拿不到全局最优，调度方案要的是可信的全局最优，不是局部可行解，所以分段线性化是工程上正确的取舍。第二个是“这套模型能不能直接执行”，答案是不能，日前调度只是计划，实际执行要配套超短期预测和滚动修正，把首时段指令下发、后续时段重优化，这一点主动说出来比被问出来强。

从那以后我每次交这类模型，不管甲方多急，都强制走一遍“约束自检脚本 + 三张图 + 成本构成表”的流程，一步不省。省下来的那点时间，往往要在答辩现场或者交付评审时十倍还回去。希望帮到你。

本文还有配套的精品资源，点击获取