企业官网建设流程全解析

先说结论：这类题目复现起来其实不难，真正劝退人的不是算法本身，而是建模细节。很多人拿到“考虑综合需求响应和碳交易机制的综合能源系统优化调度策略”这个标题，第一反应是打开 Matlab 直接写求解代码，结果半天过去连决策变量都没理清楚。我复现过几个类似的综合能源系统优化调度模型之后，最大的体会是：70% 的工作量在数学建模，20% 在数据整理，真正的求解代码反而只占 10%。这篇内容就按我实际动手的顺序，把这个题目的建模思路、Matlab 实现骨架、以及那些论文里不会写但特别容易踩的坑，完整拆一遍。

这个题目适合正在做综合能源系统方向课题的研究生、需要复现论文算法的工程师，以及想用 Matlab+Yalmip 快速搭建 MILP 调度模型的同学参考。下面直接进入正题。

1. 这题真正要复现的是什么：先看懂系统的能量流与耦合关系

1.1 综合能源系统的典型结构与设备集合

综合能源系统（Integrated Energy System, IES）之所以叫“综合”，核心在于它打通了电、热、气三种能源载体。典型的园区级系统包含这些设备：

• 热电联产机组（CHP）：最常见的是燃气轮机或内燃机，同时产电和产热，是电-热耦合的核心设备
• 燃气锅炉（GB）：补足热负荷缺口，模型比较简单，输入天然气、输出热量
• 风力发电和光伏（WT/PV）：作为清洁能源优先消纳，在调度模型中通常处理为“负负荷”或带上限的出力变量
• 电储能（ES）：电池储能，负责电负荷的削峰填谷
• 蓄热罐（HS）：热储能，给热负荷增加调节柔性
• 上级电网和天然气网：能源买入口，购电有分时电价，购气有阶梯气价

系统的能量流动可以概括为：上级电网和天然气网作为外部能源来源，经过 CHP、GB 等转换设备，配合储能和新能源，共同满足电、热、气三类负荷需求。很多复现论文还会加电转气（P2G）设备，但那是另一条技术路线，会显著增加模型复杂性，如果是第一次复现，建议先不加。

建模的第一步，是把这个系统图画成一张“能流表”：每一列是一个设备，每一行是电/热/气三种能源，交叉点标注该设备对这种能源的“产出 +1”或“消耗 -1”。这张表在后面写约束的时候非常管用，因为所有平衡约束本质上都是“某时段某能源的所有流入等于所有流出”。

1.2 综合需求响应不只是“削峰填谷”

传统电力需求响应基本只针对电负荷，比如峰谷电价引导用户错峰用电。综合需求响应（Integrated Demand Response, IDR）的关键区别在于“综合”二字：电、热、气三类负荷同时具备可调柔性，而且不同能源之间的替代关系也被利用起来。

举个最直观的例子：用户用燃气壁挂炉取暖，在电力现货价格极高的时段，可以削减一部分电负荷、同时增加一部分气负荷来补足相同的用能需求；反过来，在电价低谷时段用电采暖替代燃气供热。这种跨能源载体的负荷转移，是综合需求响应区别于普通需求响应的核心价值。

在建模层面，我通常把 IDR 拆成三类动作：

• 可削减负荷：在某一时段直接少用一部分能源，代价是按固定单价支付补偿
• 可转移负荷：把某一时段的用能需求整体平移到其他时段，需要保证一天内转移总量守恒
• 可调节负荷：负荷本身在一个区间内连续可调，比如电采暖功率允许在额定值的 40%~100% 之间波动

这三类动作分别对应不同的优化变量和约束，后文会具体展开。

1.3 碳交易机制在这类调度模型中的角色

碳交易机制从字面上看像个“环保加分项”，但在数学建模里，它直接进目标函数。这意味着它不只是让结果更好看，而是通过改变不同出力的经济成本，从机制上促进系统偏向低碳运行。

一般学术论文里对碳交易机制的建模是这套逻辑：系统根据运行情况计算实际碳排放总量，这个排放来自两部分：一是从上级电网购电对应的间接排放（用电网平均排放因子折算），二是消耗天然气对应的直接排放（按天然气排放因子折算）。系统持有一定数量的免费碳排放配额，如果实际排放低于配额，可以把富余配额拿到市场上出售获得收益；如果实际排放超过配额，超额部分需要购买配额，而且很多论文会设置阶梯式碳价——超额越多，单价越高，形成累进式惩罚。

把这个机制放入优化模型后，系统在做调度决策时会自动权衡：多买电可能便宜，但电网电的排放因子高，碳交易成本随之上升；多用气可能贵一些，但天然气的单位排放相对更可控，碳成本低。于是最优解会在“购能成本”和“碳交易成本”之间寻找平衡点。这也是为什么很多论文里只要加入碳交易机制，CHP 和燃气锅炉的出力占比就会明显变化。

2. 数学建模：从论文公式到可直接编码的约束体系

2.1 目标函数的三块成本：购能、运维、碳交易

优化调度的目标函数几乎都是总运行成本最小化，常见写法是：

min C_total = C_grid + C_gas + C_om + C_co2 + C_dr

• C_grid：向上级电网购电的成本，等于各时段购电量乘以分时电价再累加
• C_gas：购气成本，单位是元/kWh 热价或元/m³，取决于天然气的热值折算方式
• C_om：设备运行维护成本，通常是设备出力的线性函数，比如 CHP 每发一度电收 x 元维护费
• C_co2：碳交易成本或收益，如果是负值说明卖出配额获得收益
• C_dr：需求响应的补偿成本，削负荷要付补偿，转移负荷可能有激励成本

有一点必须提醒：目标函数各项的单位要完全统一。常见错误是把功率单位 MW 直接和能量单位 MWh 相加，或者购气成本用了元/m³，而设备模型里的热值折算用的是 MJ，导致最终数量级差出好几个量级。后面第 4 章会专门讲单位问题。

2.2 设备出力模型与耦合关系

设备建模是整个优化模型的地基。这里给出最常用的线性模型：

CHP 热电联产机组：

P_chp(t) = eta_chp * V_gas_chp(t) * LHV_gas % 电出力 H_chp(t) = V_gas_chp(t) * LHV_gas - P_chp(t) % 热出力由能量平衡推出

很多论文为了简化，直接使用热电比耦合约束：

H_chp(t) = K_chp * P_chp(t)

K_chp 是热电比，如果原论文没有特别说明，建议用可调区间来建模，而不是固定值：

H_chp_min <= H_chp(t) - K_chp * P_chp(t) <= H_chp_max

这样既保留了 CHP“以热定电”或“以电定热”的灵活性，又不会因为热电比参数不准确而把可行域限制死。

燃气锅炉：

H_gb(t) = eta_gb * V_gas_gb(t) * LHV_gas

燃气锅炉只产热，没有热电耦合，建模最简单。

电储能：

soc(t+1) = soc(t) + (P_ch(t) * eta_ch - P_dis(t) / eta_dis) * delta_t soc_min <= soc(t) <= soc_max P_ch(t) <= P_ch_max * u_ch(t) P_dis(t) <= P_dis_max * u_dis(t) u_ch(t) + u_dis(t) <= 1

这里需要特别注意：soc 是能量状态，而功率乘以时段长度 delta_t 才是能量变化量。如果 delta_t 是 1 小时，数值上功率和能量相等；如果后面改成 15 分钟调度（delta_t=0.25），忘记乘系数会直接导致储能状态计算错误。另外，二元变量 u_ch 和 u_dis 的互斥约束是防止储能同时充放电的必要条件，虽然惩罚项也能缓解，但最可靠的还是显式约束。

蓄热罐：建模方式与电储能完全同构，只是没有充放电互斥的约束那么严格，因为热系统惯性更大，允许一定程度的同时蓄放。

2.3 碳交易的配额分配与阶梯碳价

碳交易建模分两步：先算排放，再算成本。

实际碳排放总量的计算：

E_total = sum( P_grid(t) * EF_elec * delta_t ) + sum( V_gas(t) * LHV_gas * EF_gas )

其中 EF_elec 是电网购电的排放因子，EF_gas 是天然气燃料的排放因子。注意这里的单位换算非常讨厌：EF_elec 通常是 kgCO2/kWh，而 EF_gas 可能是 kgCO2/MJ 或者 kgCO2/m³，需要先统一到同一个能量单位体系里。

免费配额常见做法是给定一个基准排放量乘以分配比例，或者按设备额定容量折算：

E_quota = quota_ratio * E_base

超额排放量：

E_excess = E_total - E_quota

当 E_excess 为正，成本为正；当 E_excess 为负（即排放低于配额），未使用的配额可以出售，碳交易成本为负。不过很多论文为了简化只考虑超额购买的情况，复现的时候要看清原模型设计。

阶梯碳价的线性化是这一节的重点，放到第 4 章专门讲，这里先只说概念：阶梯碳价就是超额排放量越大，边际碳价越高，形成分段的线性成本函数。在 MILP 模型里，这个函数需要用分段线性化技术处理，不能直接在 Yalmip 里写 if-else。

2.4 综合需求响应的价格弹性建模

综合需求响应有两种建模路线。

路线一：价格弹性矩阵法。定义电负荷转移率与电价变化量的关系：

delta_P_e(t) / P_e_base(t) = sum_j epsilon(t,j) * delta_price(j) / price_base(j)

epsilon(t,j) 是自弹性系数（t=j）或交叉弹性系数（t≠j）。这种方法理论上很漂亮，但实际复现时有个大麻烦：弹性矩阵的参数很难从论文或实测数据里得到，很多时候论文作者自己都不清楚这些数是怎么标定出来的。我只建议在复现那些明确给出弹性矩阵数据的论文时采用这种方法。

路线二：比例边界法（推荐）。直接把需求响应能力表达为基线负荷的百分比，并给转移负荷加守恒约束：

P_load(t) = P_fix(t) + P_trans(t) - P_cut(t) sum(P_trans) = 0 % 转移负荷一天内总量守恒 -alpha * P_fix(t) <= P_trans(t) <= alpha * P_fix(t) 0 <= P_cut(t) <= beta * P_fix(t)

相比弹性矩阵法，这种方法参数更少、含义更直观——alpha 表示可转移比例，beta 表示可削减比例——而且线性约束好求解。热负荷和气负荷的建模类似，只是转移方向主要反映用能替代：比如某时段削减了电采暖负荷，同时增加了燃气供热需求。

需求响应补偿成本根据削减量和削减单价计算，转移部分如果没有明确的补偿单价约定，可以简化为不进入目标函数，只作为约束存在。

2.5 可再生能源与平衡约束

风电和光伏在调度模型里通常处理为“优先消纳”的出力变量：

0 <= P_wt(t) <= P_wt_pred(t) 0 <= P_pv(t) <= P_pv_pred(t)

如果原论文考虑“弃风弃光”，还可以在目标函数里加一项弃电惩罚，但大多数标准模型默认全额消纳，直接按预测值作为上限即可。

系统平衡约束是最后收口的三条等式：

电平衡: P_grid(t) + P_chp(t) + P_wt(t) + P_pv(t) + P_dis(t) = P_load(t) + P_ch(t) 热平衡: H_chp(t) + H_gb(t) + H_hs_dis(t) = H_load(t) + H_hs_ch(t) 气平衡: V_gas(t) + V_gas_grid(t) = V_load(t)

气网这边的平衡需要小心：有的模型把气负荷分成“燃气锅炉用气”“CHP 用气”“居民气负荷”三类，购气总量是三者之和；有的模型把 CHP 和 GB 的直接购气已经折算进设备模型。复现的时候一旦发现约束数量对不上，首先检查气平衡是否重复计算了。

3. Matlab+Yalmip 实现：变量定义、求解器配置与代码骨架

3.1 为什么选 Yalmip 而不是手写优化器

这类模型本质是一个混合整数线性规划（MILP），决策变量里有储能充放电状态这样的 0-1 变量。如果你打算用 Matlab 自带的linprog或者fmincon硬解，配置约束时非常痛苦，而且求解效率受限于算法选择。

我的经验是直接用 Yalmip 做建模层。它是一个免费的 Matlab 工具箱，核心作用是把优化变量、约束、目标函数用 Matlab 语法写出来，然后后端对接商用或开源求解器。好处很明显：

• 不需要手写约束矩阵，逻辑大幅简化
• 求解器可切换，Gurobi 不行就换 Cplex，或者切回 Matlab 自带的intlinprog
• 约束的增删改非常灵活，做情景对比时尤其方便

安装好 Yalmip 之后，确认求解器能被识别即可：

yalmiptest

3.2 决策变量的组织方式

我用 24 时段调度为例，展示一个完整的变量定义骨架。核心原则是：同一类变量定义成列向量，而不是逐时段写标量变量，这样约束和目标函数都可以向量化，代码干净且求解效率高。

T = 24; % 调度时段数，典型是 24 小时 % 设备出力变量 P_chp = sdpvar(T, 1); % CHP 电出力 H_chp = sdpvar(T, 1); % CHP 热出力 H_gb = sdpvar(T, 1); % 燃气锅炉热出力 P_wt = sdpvar(T, 1); % 风电出力 P_pv = sdpvar(T, 1); % 光伏出力 P_grid = sdpvar(T, 1); % 上级电网购电量 % 电储能 P_ch = sdpvar(T, 1); % 储能充电功率 P_dis = sdpvar(T, 1); % 储能放电功率 u_ch = binvar(T, 1); % 充电状态 u_dis = binvar(T, 1); % 放电状态 soc = sdpvar(T+1, 1); % SOC，T+1 是为了包含初始时刻和结束时刻 % 需求响应变量 P_trans = sdpvar(T, 1); % 可转移电负荷（正为转入负荷） P_cut = sdpvar(T, 1); % 可削减电负荷 % 碳交易相关 E_total = sdpvar(1, 1); % 总碳排放量 C_co2 = sdpvar(1, 1); % 碳交易成本

binvar对应 0-1 变量，sdpvar对应连续变量。储能变量的 T+1 维度处理是很多新手容易忽略的——SOC 的状态转移约束会涉及 t=1 和 t=T+1 两个边界时刻。

3.3 目标函数与约束的代码化

目标函数建议按成本项分块计算，最后求和。这样后续调试时可以单独查看每一项的值，方便定位异常。

% 分时电价、气价、运维系数、碳排放参数（示例数据） price_e = [0.35; 0.5; ...; 0.3]; % 24x1，元/kWh price_g = 0.32; % 元/kWh（天然气热值折算后） delta_t = 1; % 时段长度 1h % 目标函数分项 C_grid = price_e' * P_grid * delta_t; % 购电成本 C_gas = price_g * V_gas_total * delta_t; % 购气成本 C_om = k_chp' * (P_chp + H_chp) + k_gb' * H_gb; % 运维成本 C_dr = dr_price' * P_cut * delta_t; % 需求响应补偿 objective = C_grid + C_gas + C_om + C_co2 + C_dr;

约束部分，我用一个循环来写储能的状态转移约束，这是最直观的写法：

Constraints = []; % 储能SOC递推约束 for t = 1:T Constraints = [Constraints, soc(t+1) == soc(t) + ... (P_ch(t) * eta_ch - P_dis(t) / eta_dis) * delta_t]; end % 初始和最终SOC约束 Constraints = [Constraints, soc(1) == soc_0]; Constraints = [Constraints, soc(T+1) == soc_0]; % 周期性运行要求 % 储能充放电互斥 Constraints = [Constraints, P_ch <= P_ch_max .* u_ch]; Constraints = [Constraints, P_dis <= P_dis_max .* u_dis]; Constraints = [Constraints, u_ch + u_dis <= 1];

注意等式约束用的是==，Yalmip 对浮点误差很敏感。如果你的参数导致矩阵病态，经常会出现“不可行”或者求解结果明显违背物理规律，优先检查是不是系数差了好几个数量级。

其他约束，比如功率平衡、设备上下限、需求响应变量边界，写法上完全一致。把所有的Constraints累积起来，最后调用求解器：

opts = sdpsettings('solver', 'gurobi', 'verbose', 2); ops = optimize(Constraints, objective, opts);

求解完成之后，取值用value()：

P_chp_opt = value(P_chp); P_grid_opt = value(P_grid); soc_opt = value(soc);

3.4 求解器选型与参数设置

求解器选择直接影响能不能解出来、解多快。我按优先级排序：

对于 24 时段的园区级 IES 调度模型，决策变量数通常不到 500 个（含 0-1 变量），这属于小规模 MILP，intlinprog也能在几秒到几十秒内求解。但如果你想做 96 时段（15 分钟）或者含多场景随机优化的扩展，建议还是装 Gurobi。sdpsettings里比较实用的几个参数也提一下：

opts = sdpsettings('solver', 'gurobi', ... 'verbose', 1, ... % 打印迭代日志 'gurobi.MIPGap', 0.0001, ... % MIP相对间隙，越小越精确 'gurobi.TimeLimit', 300); % 300秒限制，防止卡死

4. 复现路上最容易翻车的五个细节

4.1 单位与数量级错位

这是我在帮别人排查代码时遇到最多的问题。论文里经常混用这几类单位：功率（kW/MW）、能量（kWh/MWh）、热值（MJ/m³、kWh/m³）、排放因子（kgCO2/kWh、kgCO2/MJ）。

一个典型事故现场：购气价格是 2.8 元/m³，而设备出力模型用的是 kW 和 kWh，天然气的低位热值按 9.97 kWh/m³ 算。如果你直接把购气量乘气价，但设备那边把 m³ 当成了 kWh，就会把燃气能量算大或算小。调度结果会出现“CHP 疯狂发电、电储能完全不充电”这种明显不合理的现象。

我的建议是：在写代码之前，把所有参数统一折算到一个单位体系里。我习惯统一使用 MW 和 MWh，天然气能量统一用 MWh，碳排放统一用 tCO2（吨）。在代码里专门建一个参数文件，把原始参数和统一后的参数对照写清楚：

LHV_gas = 9.97; % kWh/m3 LHV_gas_MWh = LHV_gas / 1000; % MWh/m3 gas_price_per_kwh = 2.8 / LHV_gas; % 元/kWh 热价 gas_price_per_mwh = gas_price_per_kwh * 1000; % 元/MWh 热价

最后检查目标函数的量纲：所有成本项加起来应该单位一致，都是元，这样目标函数才有意义。

4.2 储能 SOC 的初值/终值约束

如果没有给 SOC 加初始值和终值约束，会出现一个很搞笑的优化结果：储能系统在最后一个时段把电全部放光，SOC 直接降到下限，因为“偷掉”这个能量可以减少购能成本。这在单日调度里是可行解，但不是实际可运行的策略。

解决办法是在模型里要求系统按日循环运行：

Constraints = [Constraints, soc(1) == 0.2]; % 初始SOC Constraints = [Constraints, soc(T+1) == soc(1)]; % 周期回归

如果论文中没有明确给出初始 SOC，我一般取 SOC 上下限的中值，然后强制终值等于初值。这样储能的作用是削峰填谷，而不是额外提供“免费能量”。

4.3 阶梯碳价的分段线性化

很多论文里的碳交易机制是阶梯累进碳价：超额排放量在第一个区间内的碳价是 c1，超过一定量后进入第二个区间，碳价跳到 c2，以此类推。这种分段函数在 Yalmip 里不能用 if-else，需要显式引入 0-1 变量来做区间指示。

假设超额排放量 x 有三个区间[0,a]、[a,b]、[b,+inf)，对应累进碳价 c1<c2<c3，可以这样线性化：

a = 100; b = 200; % 区间边界 c1 = 60; c2 = 90; c3 = 120; % 单位碳价（元/吨） M = 1e4; % 大M常数 u1 = binvar(1,1); u2 = binvar(1,1); u3 = binvar(1,1); x1 = sdpvar(1,1); x2 = sdpvar(1,1); x3 = sdpvar(1,1); Constraints = [Constraints, x == x1 + x2 + x3]; Constraints = [Constraints, 0 <= x1 <= a * u1]; Constraints = [Constraints, 0 <= x2 <= (b - a) * u2]; Constraints = [Constraints, 0 <= x3 <= M * u3]; Constraints = [Constraints, u1 + u2 + u3 == 1]; C_co2 = c1 * x1 + c2 * x2 + c3 * x3;

这样碳交易成本就变成了一个线性表达式，可以安全放进 MILP。大M的值不要取得过大，取略大于该变量物理上限的数量级即可，太大会导致数值求解困难。

4.4 需求响应变量的可行域

综合需求响应建模时，转移负荷最容易出问题。忘记写sum(P_trans) == 0这个守恒约束，优化器会直接在低电价时段疯狂“转入”负荷、在高电价时段“转出”负荷，人为拼出无穷的能量来源，最终得到一个失真且毫无意义的最优解。

另外，可削减负荷的上下限也要和基线负荷挂钩。如果没有显式约束，优化器会把削减量开到无穷大来降低成本。正确写法是：

Constraints = [Constraints, sum(P_trans) == 0]; Constraints = [Constraints, - alpha .* P_fix <= P_trans <= alpha .* P_fix]; Constraints = [Constraints, 0 <= P_cut <= beta .* P_fix];

这里的 alpha 和 beta 是论文给的负荷可调比例。如果论文没给，常见范围是 0.1 到 0.2，我一般先按 0.1 跑通模型，再扫参看结果趋势。

4.5 情景对比实验的设计

复现论文还有一个隐藏任务：复现论文里的实验结果图，尤其是那几个对比情景的图表。如果只看最终最优解，不做情景对比，很难说明需求响应和碳交易机制各自的作用。

标准的消融实验设计是四个情景：

每个情景除了这两个开关外，其他所有参数（设备参数、负荷数据、分时电价、气价、风电光伏预测）必须完全一致。我在代码里会用一个标志位控制模型生成逻辑：

flag_idr = 1; % 是否考虑综合需求响应 flag_co2 = 1; % 是否考虑碳交易机制 if ~flag_idr % 不添加负荷转移和削减的可调变量，直接把负荷当成固定值 Constraints = [Constraints, P_load == P_fix]; else Constraints = [Constraints, P_load == P_fix + P_trans - P_cut]; % ... 其他需求响应约束 end if flag_co2 % 目标函数里加入 C_co2 objective = objective + C_co2; end

这样跑四次，收集四组结果，后面画对比图就非常顺手。

5. 结果可视化与灵敏度分析该怎么写

5.1 设备出力与功率平衡图

复现类工作最直观的成果就是调度计划曲线，很多人一到画图环节就把解锁代码写成一堆散点，看起来乱。

我通常用stairs画出阶梯状曲线，因为 24 时段离散调度的本质就是阶梯式决策，stairs比plot更真实。核心设备出力画在同一张图上：

figure; stairs(1:T, P_grid_opt, 'LineWidth', 1.5); hold on; stairs(1:T, P_chp_opt, 'LineWidth', 1.5); stairs(1:T, P_wt_opt + P_pv_opt, 'LineWidth', 1.5); stairs(1:T, P_dis_opt - P_ch_opt, 'LineWidth', 1.5); legend('购电','CHP','新能源','储能净放电','Location','best'); xlabel('时段/h'); ylabel('功率/MW'); grid on;

再画一张功率平衡验证图，把电负荷曲线和所有电源出力叠加，确认每个时段加总都能对上。这一步虽然是验证步骤，但很多人忽略，导致论文里出现“平衡约束没满足”的低级错误。

5.2 需求响应前后的负荷曲线对比

需求响应的效果一般用原始负荷曲线和优化后的实际负荷曲线对比来呈现，可以直观看到削峰填谷的效果：

figure; stairs(1:T, P_fix, '--', 'LineWidth', 1.5); hold on; stairs(1:T, value(P_load), 'LineWidth', 1.5); legend('原始负荷','需求响应后负荷','Location','best'); xlabel('时段/h'); ylabel('负荷/MW'); grid on;

同时可以把四个情景的总成本、碳排放量、新能源消纳量整理成一个结果汇总表，这是论文结论部分的标准呈现方式。用表格输出时，四舍五入的精度要注意，成本列保留两位小数，碳排放列保留三位比较合适。

5.3 碳交易价格与配额等参数敏感性分析

灵敏度分析是这个题目里最容易出彩的部分。很多论文在算例分析部分都有这么一张图：横轴是碳交易价格，纵轴是系统总成本和碳排放量，观察两条曲线如何相交或变化。

我自己复现时会把核心参数封装成函数，然后写一个循环扫描碳价：

carbon_prices = 20:10:200; total_cost_list = zeros(size(carbon_prices)); co2_list = zeros(size(carbon_prices)); for i = 1:length(carbon_prices) carbon_price = carbon_prices(i); [total_cost_list(i), co2_list(i)] = run_ies_model(carbon_price); end figure; yyaxis left; plot(carbon_prices, total_cost_list, 'o-', 'LineWidth', 1.5); ylabel('系统总成本/元'); yyaxis right; plot(carbon_prices, co2_list, 's-', 'LineWidth', 1.5); ylabel('碳排放量/吨'); xlabel('碳价/(元/吨)'); grid on;

这种参数扫描做起来很快，而且能讲清楚一个关键结论：随着碳价升高，系统的碳排放量下降，但总成本不一定单调上升。因为需求响应和碳交易相互配合时，系统会通过调整负荷曲线和能源结构，把高碳出力的时段避开，部分成本被负荷转移带来的购能成本下降抵消。

如果在 5.3 之后要自然收尾，我一般会加一小段个人操作体会。比如，我在跑完碳价敏感性之后发现，如果只开需求响应不开碳交易，系统的总成本下降有限，但负荷曲线明显变平；如果只开碳交易不开需求响应，碳排放降得很明显但成本上升更快。两个机制同时开启时，才存在一个成本可以接受同时排放也有下降的甜点区间。这个观察不一定所有论文都明确写出来，但对理解模型机理非常关键。在参数扫描时多花十几分钟把这类交叉现象找出来，往往比复制论文的图更有价值。

最后再分享一个小技巧：这类模型调试时如果遇到“infeasible problem”，不要急着去翻约束逐条看，先用 Yalmip 的check(Constraints)找出哪些约束的值是 NaN 或 Inf，然后专门盯着那几条约束的变量量纲查。绝大多数不可行问题，根源都是单位不统一导致某条平衡约束的系数差了 10 的倍数。把这一关过了，剩下的求解和出图通常就顺畅很多。