★ FEATURED ARTICLE 做综合能源系统优化运行的人应该都遇到过这种局面电价在高峰时段飙升燃气轮机满发还是不够用用户侧的负荷峰谷差又大得离谱等你把机组出力优化完调度一看又冒出来一句话——碳排放超标了。这个时候你才会意识到只盯着供能侧做经济调度路子是走窄了。负荷侧明明有弹性空间碳排指标也早就变成了真金白银的成本项这两块不参与优化整个模型就是瘸腿的。这篇要聊的就是一套完整的MATLAB实现方案在碳交易机制下把需求响应和柔性负荷纳入综合能源系统的优化运行调度。所谓碳交易机制下考虑需求响应的综合能源系统优化运行翻译成大白话就是——供能侧有燃气轮机、燃气锅炉、电储能、热储能这些设备负荷侧有可转移、可削减、可平移的柔性负荷系统在满足电热负荷需求的前提下把购电、购气、运维、碳交易、需求响应补偿这几笔账放在一个目标函数里统一优化最终给出各设备逐时出力曲线和总成本最低的运行方案。这套东西非常实用适合做能源系统优化、微电网调度、碳减排路径分析的研究生和工程师参考也可以直接作为论文或课题的基线模型。下面我会按照模型为什么这样设计→供能侧数学建模→碳交易与柔性负荷的约束写法→MATLAB代码链路→调参踩坑这条线把整个逻辑完整拆开讲。1. 为什么要把碳交易和需求响应放进同一个优化模型1.1 综合能源系统的电热耦合账不能分开算先说清楚综合能源系统到底综合在哪里。常规微电网只算电冷热负荷往往单独用锅炉或者空调扛。但真正的综合能源系统典型构成是燃气轮机CHP、燃气锅炉、电储能、热储能这几样电和热通过CHP机组天然耦合在一起。燃气轮机在生产电能的同时产生高温烟气余热可以用来供热供电和供热本质上是同一份天然气消耗的两种产出。这意味着什么意味着你优化供电时必须同时决定供热方式。比如晚间电价便宜我可以让CHP少发电多用电网购电热负荷缺口用燃气锅炉补但如果燃气价格便宜而电价贵最优策略可能反过来——CHP多发电余热满足一部分热负荷锅炉少烧气。这种决策在数学上就是一个耦合优化问题单独算电或者单独算热都会得到次优甚至不可行的调度方案。我在实际搭建模型时会把电功率平衡和热功率平衡作为两个最核心的等式约束其他所有设备、储能、需求响应都挂在它们下面。1.2 碳交易机制碳排放变成了有价格的生产要素碳交易的本质是给二氧化碳排放量定了一个市场价。系统每消耗一份天然气、每从电网买一度电如果电力来自高碳电源都对应着一定碳排放。碳交易机制通常先给系统分配一个免费配额实际排放低于配额富余部分可以卖钱超过配额就要购买碳排放权而且超过越多单位碳价越贵——这就是常说的阶梯碳价。我给模型加入碳交易机制后最直观的体会是燃气轮机不再多发电就是好因为发电越多天然气消耗越大碳排放成本越高。原本经济调度中最优的电转热储热与储电配合策略在碳价压力下可能完全逆转。碳交易本质上是在供能侧加了一个价格杠杆而这个价格杠杆会直接影响设备出力的排序。1.3 需求响应与柔性负荷负荷侧也是可调资源很多做调度的人有一个思维惯性负荷曲线是给定的是硬约束优化任务只是让发电/供能去追负荷。但真实系统的负荷侧是有弹性的。工业用户可以把部分工序从高峰时段挪到低谷商业楼宇的空调可以在短时间内降功率运行而用户几乎感知不到居民侧的洗衣机和电动车充电时间也可以平移。这些统称为柔性负荷。把它们纳入优化模型后负荷曲线就不再是外生给定的而是部分内生的优化变量。系统在电价高的时段请求用户少用一点电在电价低的时段鼓励用户多用一点电整个系统的负荷峰谷差可以大幅改善。当然让用户配合调整不是免费的要么按激励价格给予补偿要么通过分时电价引导——这笔钱要进目标函数。这正好与碳交易形成呼应碳交易影响供能侧成本需求响应重塑负荷侧曲线两者共同作用才是一个完整的、含柔性负荷的综合能源系统优化模型。2. 供能侧建模设备、储能与电热平衡约束的完整展开2.1 燃气轮机和燃气锅炉耗能转换设备的数学描述设备建模是整个优化模型的地基。我用的是典型的能源集线器结构核心设备是CHP燃气轮机和燃气锅炉。CHP机组有两个输出——电功率 (P_{chp}(t)) 和热功率 (Q_{chp}(t))消耗天然气 (F_{chp}(t))。如果不引入复杂的变工况效率曲线工程中常用分段常数效率来近似[ P_{chp}(t) \eta_e \cdot F_{chp}(t), \quad Q_{chp}(t) \eta_h \cdot F_{chp}(t) ]实际建模时我会把两个式子合并成热电比约束比如 (Q_{chp}(t) k_{chp} \cdot P_{chp}(t))也就是热出力与电出力固定比例。这样约束少了模型也稳。如果追求精度也可以让热电比在一个区间内变化但代价是引入更多非线性项求解难度会明显上升。燃气锅炉相对简单就是一个天然气进、热水出的转换器[ 0 \le Q_{gb}(t) \le Q_{gb,max}, \quad C_{fuel,gb}(t) \frac{Q_{gb}(t)}{\eta_{gb}} \cdot c_{gas} ]设备建模还需要上下限约束和爬坡约束。爬坡约束尤其容易被忽略但它对真实系统非常重要——燃气轮机出力不可能一分钟内从10%跳到100%。我通常写成[ -P_{chp}^{ramp} \le P_{chp}(t) - P_{chp}(t-1) \le P_{chp}^{ramp} ]注意这里的索引是从第2个时段开始到T循环别写错后面调试篇我会专门讲这种坑。2.2 电储能与热储能动态约束和周期性问题储能是综合能源系统里最灵活的缓冲器。电储能的状态方程是[ E_{es}(t1) E_{es}(t) \cdot (1 - \sigma_{es}) P_{ch,es}(t) \cdot \eta_{ch} - \frac{P_{dis,es}(t)}{\eta_{dis}} ]同一时段内不能同时充电和放电这个约束在工程上通常用二进制变量实现[ P_{ch,es}(t) \le M \cdot u_{ch}(t), \quad P_{dis,es}(t) \le M \cdot u_{dis}(t), \quad u_{ch}(t) u_{dis}(t) \le 1 ]M要取得足够大但又不要过大否则会拖慢求解器速度。我自己一般取设备额定功率的1.2倍。热储能类似状态变量是蓄热罐的热量单位是kWh或MWh公式结构完全一样。储能建模有一个细节运行周期结束后要不要恢复到初始状态。如果是研究典型日调度我通常加上 (E_{es}(1) E_{es}(T1))让储能日循环闭合避免模型通过偷用初始储能来压低成本——这个问题不做会吃大亏。2.3 电、热功率平衡整个模型的骨架不管前面加了多少设备最终都要汇聚在功率平衡约束上。电平衡方程我写成[ P_{buy}(t) P_{chp,e}(t) P_{dis,es}(t) P_{load}(t) P_{ch,es}(t) P_{shift,in}(t) - P_{shift,out}(t) - P_{cut}(t) ]注意柔性负荷进入平衡约束的方式是把负荷侧从固定值改成可变量——可转移负荷转出的部分相当于负荷减小转入的部分相当于负荷增加可削减负荷被削减的部分直接视为负荷下降。把这些都放在等式右边而不是在左边加负的发电逻辑上更清晰也方便后续输出结果时拆分统计。热平衡约束相对简单[ Q_{chp,h}(t) Q_{gb}(t) Q_{dis,hs}(t) Q_{load}(t) Q_{ch,hs}(t) ]这两条约束是整个模型运行逻辑的骨架所有设备的出力决策最终都受这两条等式约束牵制。我调试模型时第一步永远是把各种成本置零、只留设备上下限然后检查这两条平衡约束是否在每一个时段都严格满足。3. 碳配额、碳价阶梯与柔性负荷的数学化表达3.1 实际碳排放与免费配额的计算口径碳交易建模的关键在于两套口径实际碳排放量和免费配额。实际碳排放量我这里按天然气消耗量换算[ E_{CO2}(t) k_{CO2} \cdot \left( F_{chp}(t) F_{gb}(t) \right) ]其中 (k_{CO2}) 是天然气燃烧的碳排放因子单位是 tCO2/MWh。不同文献取法略有差异我常用的近似值是0.20左右具体数值按你采用的燃料热值和排放因子折标。从电网购电是否算入碳排放取决于你的模型边界——如果研究的是园区系统外购电力的碳排通常按电网平均排放因子单独计如果电网侧碳责任不在本系统内可以不含这项。这一点最好在论文里写明口径否则评审会挑刺。免费配额我用基准线法[ M_{quota} \delta_e \cdot P_{chp,e}(t) \delta_h \cdot \left( Q_{chp,h}(t) Q_{gb}(t) \right) ](\delta_e) 和 (\delta_h) 分别是供电和供热的配额基准值。配额的逻辑是按实际供能量乘以行业基准线折算也就是说你提供了多少能量服务就免费获得多少排放空间。3.2 阶梯碳价的分段线性化用二进制变量实现实际碳排放减去免费配额就是需要购买碳权的净排放量。为了体现超排越多罚得越狠的政策导向碳价通常设计成阶梯函数。我常用的分段方式是若 (E - M \le 0)配额富余富余部分按出售价 (c_{sell}) 获得收益若 (0 \le E - M \le d_1)超排部分按基准碳价 (c_0) 购买若 (d_1 E - M \le d_1 d_2)超排部分单价从 (c_0) 涨到 (c_1)若超过 (d_1 d_2)继续按更高价 (c_2) 购买。这种分段函数直接用 if-else 在MATLAB里表达很简单但Yalmip无法直接处理需要用二进制变量把区间选出来。我的做法是引入非负连续变量 (\Delta_1, \Delta_2, \Delta_3)分别代表落在第一、二、三超排区间的排放量然后加约束( \Delta_1 \Delta_2 \Delta_3 \max(E - M, 0) )( 0 \le \Delta_1 \le d_1 )( 0 \le \Delta_2 \le d_2 \cdot z_1 )( 0 \le \Delta_3 \le d_3 \cdot z_2 )其中 (z_1, z_2) 是二进制变量含义是是否进入了第二/第三区间。只有前面区间用满了后面区间才能取非零值——这个前后联动需要通过约束保证比如 (\Delta_2 \ge d_2 \cdot z_2) 这样的逻辑。我最开始写这里时吃过亏后面调试章节会具体展开。3.3 三类柔性负荷的详细约束写法柔性负荷我分成三类来建模这是整个模型中需求响应的核心细节。可转移负荷典型例子是电动汽车充电、洗衣机、蓄冷蓄热设备。特征是在一定调度周期内总用电量不变但用电时段可以整体或局部平移。约束写法[ \sum_{t1}^{T} P_{shift}(t) 0, \quad -\overline{P_{shift}} \le P_{shift}(t) \le \overline{P_{shift}} ]这里 (P_{shift}(t)) 我用正负号表示转入/转出功率。总转移量为0保证只是挪时间不是凭空增减用电。可削减负荷典型例子是空调温度调节、照明减半。特征是可以中断或降低一部分功率用户会有一定舒适度损失系统需要给予补偿。约束写法[ 0 \le P_{cut}(t) \le \overline{P_{cut}}(t), \quad \sum_{t1}^{T} P_{cut}(t) \le \overline{P_{cut,total}} ]第二种约束是总削减量比例上限防止模型为了让成本最低而不停削负荷把用户体验完全牺牲掉。可平移负荷典型例子是工厂的固定工序流程。这类负荷整体在时间轴上滑动比如一条生产线持续运行4小时可以选择从早上8点开始也可以从上午10点开始。建模比前两类复杂需要引入整数变量表示启动时段[ \sum_{\tau0}^{L-1} P_{move}(t\tau) P_{base,profile} \cdot \sum_{s \in S} \delta_{s, t} ]这里工程上用二进制变量 (\delta_s) 表示在时段 s 启动该负荷然后整条负荷曲线随之平移。如果觉得这个建模太繁琐初期可以先用可转移负荷替代效果类似但实现简单很多。三种类型负荷的激励成本不一样可转移负荷通常走分时电价引导可削减负荷需要按单位功率给补偿可平移负荷可能是签订合同后的固定服务费。这个差异要体现在目标函数里下一节会展开。3.4 目标函数四类成本怎么组装整个优化模型的目标函数我写成[ \min ; C_{total} C_{grid} C_{gas} C_{om} C_{carbon} C_{dr} ](C_{grid})从电网买的电费按分时购电价与购电功率乘积在T个时段上求和(C_{gas})天然气费用CHP和锅炉的耗气量之和乘以气价(C_{om})设备运行维护成本一般按出力线性计费每类设备一个单位运维系数(C_{carbon})阶梯碳交易成本按上一节的分段函数计算(C_{dr})需求响应成本主要是可削减负荷的补偿费用再加上可平移负荷的合同费用。这个目标函数组装起来非常直观Yalmip里就是把每一部分的表达式加起来。但有个细节值得注意碳交易成本如果是阶梯函数它在目标函数里会引入非凸的分段线性项。如果分段数不多直接用二进制变量加辅助变量的方法处理即可如果分段数很多建议用Yalmip内置的implies逻辑否则约束数量会爆炸。4. MATLAB代码实现从参数准备到结果输出4.1 代码结构设计把模型拆成清晰模块我通常把整个项目的MATLAB代码拆成五个文件main_IES_DR_CET.m % 主脚本负责整体流程 config_params.m % 参数定义设备容量、效率、价格、碳排放因子、柔性负荷参数 build_decision_vars.m % 决策变量定义连续变量 二进制变量 build_constraints.m % 约束组装设备约束、储能耗约束、平衡约束、碳阶梯、柔性负荷 build_objective.m % 目标函数组装 solve_and_plot.m % 求解、结果提取、绘图与表格输出main_IES_DR_CET.m的核心代码逻辑如下config_params; % 第1步加载所有参数 build_decision_vars; % 第2步创建sdpvar变量 build_constraints; % 第3步逐段添加约束 build_objective; % 第4步定义目标函数 ops sdpsettings(solver, cplex, verbose, 2); result optimize(Constraints, Objective, ops); solve_and_plot; % 第5步提取变量并画图这种结构对调试非常友好。比如想单独验证没有碳交易的基准场景只要在build_constraints里注释掉碳阶梯约束、在build_objective里去掉碳成本项即可不用动其他文件。4.2 用Yalmip声明变量和约束的示例决策变量我用sdpsvar声明核心代码长这样T 24; % 24时段调度 % CHP机组 P_chp sdpvar(1, T, full); Q_chp sdpvar(1, T, full); % 燃气锅炉 Q_gb sdpvar(1, T, full); % 电储能 E_es sdpvar(1, T1, full); P_ch_es sdpvar(1, T, full); P_dis_es sdpvar(1, T, full); u_es binvar(1, T, full); % 1表示充电 % 柔性负荷 P_shift sdpvar(1, T, full); % 可转移净转入功率正为转入 P_cut sdpvar(1, T, full); % 可削减功率注意电储能状态变量定义成T1维因为递推关系需要用到末时段状态。约束组装部分第一个核心约束是爬坡for t 2:T Constraints [Constraints, -P_chp_ramp P_chp(t) - P_chp(t-1) P_chp_ramp]; end电平衡约束Constraints [Constraints, ... P_buy P_chp P_dis_es P_load_base P_ch_es P_shift - P_cut];这行写得很紧凑但注意MATLAB向量化加法要求所有变量都是一样维度P_load_base必须是你输入的24维负荷向量P_shift正负号已经表示了转出转入方向。如果出现维度报错十有八九是某一行变量维度不统一。碳交易约束的实现可以先算辅助变量E_carbon k_co2 * (F_chp F_gb); % 实际碳排放 M_quota delta_e * P_chp delta_h * (Q_chp Q_gb); % 免费配额 E_net E_carbon - M_quota; % 净排放可能为负然后针对净排放的正负和阶梯区间加判断。Yalmip里直接用max(E_net, 0)在大多数情况下也可以处理但如果要求严格的分段线性建议手动引入二进制变量。我在第3.2节已经给了约束框架这里不再重复。4.3 设置求解器、结果输出与场景对比求解设置这一步最简单但也很容易出问题ops sdpsettings(solver, cplex, verbose, 1, showprogress, 1); diagnostic optimize(Constraints, Objective, ops); if diagnostic.problem 0 disp(求解成功); P_chp_opt value(P_chp); P_shift_opt value(P_shift); P_cut_opt value(P_cut); ... else disp(求解失败); disp(diagnostic.info); end拿到最优解以后我强烈建议把所有关键指标都算出来存成表格。下面是我常用的三种对比场景场景说明关注指标基准场景不含DR、不含碳交易总成本、碳排放仅碳交易含阶梯碳价不含DR碳成本、机组出力变化碳交易DR含碳交易和全部柔性负荷总成本、峰谷差、削减量有了这三个场景的数据你才能说明碳交易加需求响应这个组合到底带来了多少收益。我实际跑出来的典型结果是加入碳交易后CHP机组在高碳价时段出力会降低燃气锅炉和电储能承担更多调节压力再加入需求响应后高峰时段负荷被转移和削减购电成本明显下降碳成本也随之下降总成本比基准场景低不少。这些结论在论文里就是很好的经济性与环保性分析。5. 调试与踩坑跑通模型时最容易翻车的几个地方5.1 储能SOC的维度错位问题电储能状态E_es定义成T1维但控制变量P_ch_es和P_dis_es是T维。我一开始写递推约束时循环了t 1:T结果MATLAB报索引越界或者更隐蔽的是E_es(t1)被当成未知变量但从未出现在最后的输出里导致模型变量数比约束多出现冗余自由度。解决办法很简单一定要检查约束中所有变量的维度并且最后求解完用value(E_es(1))对比初始SOC看是不是你设定的E_init。如果初始值对不上大概率是递推关系写错了。5.2 阶梯碳价区间没用满引发的约束错乱这是碳交易建模最容易出的问题。我刚写阶梯碳价时直接约束 (0 \le \Delta_1 \le d_1)、(0 \le \Delta_2 \le d_2) 就完事了结果优化器永远只让超排量落在第一区间因为没有约束强制进入第二区间前必须先用满第一区间。这会让阶梯惩罚彻底失效碳价形同虚设。加一个逻辑约束就好了Constraints [Constraints, implies(z_2, Delta_1 d_1)];含义是如果启用了第二区间那么第一区间必须正好用满。类似的第三区间启用时第二区间也必须是满的。加上这一组约束之后阶梯碳价才是真正的阶梯。5.3 求解超时与不可行解的定位方法模型规模不大时Cplex秒解但一旦加入多个二进制变量储能不能同时充放、可平移负荷启动标志、碳区间选择求解时间可能显著上升。我的经验是先跑一个设备无爬坡、无储能状态、无碳阶梯、无DR的极简版本验证线性模型本身没问题逐步加约束每加一类就检查结果是否合理如果出现不可行解用diagnostic.info看具体报错类别然后逐个注释疑似有问题的约束来定位。另外数据量纲不统一是常见的隐性问题。电价是元/kWh气价是元/m³碳排放因子是tCO2/MWh——三套单位混在一起目标函数里两项的成本量级可能差了1000倍。建议全部统一到MWh、元/MWh量纲否则目标函数里数值大的项会压制数值小的项优化器看不见小量级成本。5.4 给新手的实操建议从哪里开始改代码如果你拿到一份类似代码想跑通我建议按照这个顺序来先把config_params.m里的参数全部改成你自己的算例数据尤其是负荷曲线和电价曲线跑一个什么都不加的基准场景确认平衡约束没爆加储能检验SOC曲线是连续合理的加碳交易检验不同碳价下CHP出力是否符合直觉最后加需求响应重点看可转移负荷的24小时净转移功率之和是否为0可削减负荷是否超过了设定的总上限。我实际调试中发现绝大多数人卡在最后一步——柔性负荷约束加进去以后平衡方程直接不可行解决方法是检查可转移功率的最大值是否给了合理上界。如果你让P_shift上限设得过大模型会试图把所有负荷都挪走而平衡方程的物理意义就崩了。写在最后的一点心得这套系统我断断续续改过很多版最深的体会是模型复杂度一定要和你的目标匹配。如果是写论文碳交易分层三档、柔性负荷三类全上完全可以撑起一篇核心章节如果是项目预研其实可转移负荷加一档阶梯碳价就够用跑得又快又稳。另外参数敏感性分析非常值得做——把碳价基准值从40元/吨扫到100元/吨把需求响应补偿倍数从0.8扫到2.0看系统总成本和碳排放的变化趋势这类图放在报告里很有说服力也能帮你验证模型行为是否符合真实系统的经济学逻辑。跑模型的时候记得随手保存每个场景的workspace变量后面画对比图和写分析段落会省很多时间。
阅读完成 · 觉得有帮助？
