考虑阶梯式碳交易机制与电制气的综合能源系统热电优化

参考文献：《考虑阶梯式碳交易机制与电制气的综合能源系统热电优化》
摘要:“双碳”背景下,为提高能源利用率,优化设备的运行灵活性,进一步降低综合能源系统(IES)的碳排放水平,提出一种IES低碳经济运行策略。首先考虑IES参与到碳交易市场,引入阶梯式碳交易机制引导IES控制碳排放:接着细化电转气(P2C)的两阶段运行过程,引入电解槽、甲烷反应器、氢燃料电池(HFC)替换传统的P2C,研究氢能的多方面效益:最后提出热电比可调的热电联产、HFC运行策略,进一步提高IES的低碳性与经济性。基于此,构建以购能成本、碳排放成本、弃风成本最小的低碳经济运行目标,将原问题转化为混合些教线性问题,运用CPLEX商业求解器进行求解,通过设置多个运行情景,对比验证了所提策略的有效性。
关键词:氢能;阶梯式碳交易机制;热电比可调;综合能源系统;低碳经济

1 引言

社会经济发展给环境带来了巨大的负担口,温室气体的大量排放进一步导致全球气候变暖口。在此大背景下,我国表示力争在 2030 年前实现碳达峰.2060年前实现碳中和国。电力行业作为能源消耗的主体,其碳排放量在碳排放总量中占了很大的比重实现低碳电力将有望加速实现碳减排的目标。综合能源系统IES(Integrated Energy System)内部耦合了多种能源进行联合供应,能满足终端多能负荷需求,进一步优化了多能系统的低碳经济性。

综合能源系统是指将多种能源资源（如化石能源、可再生能源等）进行整合和优化配置，以满足能源需求、提高能源利用效率和降低环境影响的系统。综合能源系统的核心理念是将各种能源资源进行协调和协同，实现能源的高效利用和综合优化，从而实现能源供应的可持续性和稳定性。
综合能源系统的特点包括多元化能源资源、高效能源利用、灵活性和可靠性。通过整合多种能源资源，综合能源系统可以更好地应对能源供应的不稳定性和波动性，提高能源供应的可靠性和稳定性。同时，综合能源系统可以实现能源资源的互补利用，提高能源利用效率，降低能源消耗和碳排放。
综合能源系统包括多种技术和设施，如多能互补系统、能源存储技术、智能能源管理系统等。多能互补系统是指将不同能源资源（如太阳能、风能、生物质能等）进行整合和协同，以满足能源需求。能源存储技术可以将能源在时间和空间上进行储存和调度，提高能源供应的灵活性和可靠性。智能能源管理系统可以通过数据分析和优化算法，实现对能源系统的智能监控和调度，提高能源利用效率和系统运行的稳定性。
综合能源系统在能源转型和可持续发展方面具有重要意义。通过整合和优化各种能源资源，综合能源系统可以实现能源供应的多元化和灵活性，提高能源利用效率，降低碳排放和环境影响。同时，综合能源系统也可以促进能源技术的创新和应用，推动能源产业的转型升级，实现能源供应的可持续性和稳定性。
综合能源系统是一种创新的能源供应模式，具有提高能源利用效率、降低环境影响和推动能源转型的重要作用。通过不断完善和推广综合能源系统，我们有望实现能源供应的可持续性和稳定性，推动能源产业的绿色发展，建设更加清洁和美丽的环境。
阶梯碳交易是一种碳排放交易制度，旨在通过设定不同阶梯的碳排放标准来鼓励减少碳排放和提高能源效率。这种制度通过对不同行业和企业设置不同的碳排放标准，激励企业采取更加环保的生产方式，减少碳排放，从而实现减排目标。在阶梯碳交易制度中，碳排放标准通常以每单位产出的碳排放量来衡量，不同行业和企业根据其产业特点和规模设置不同的碳排放标准。一般来说，排放量较高的企业需要支付更多的碳排放费用，而排放量较低的企业则可以获得碳排放配额的奖励。
阶梯碳交易制度的实施可以有效地激励企业减少碳排放，促进绿色发展。首先，通过设置不同的碳排放标准，可以引导企业采取更加环保的生产方式，减少碳排放。其次，通过对碳排放费用的收取和奖励机制的建立，可以激励企业主动减排，提高能源利用效率，降低碳排放成本。另外，阶梯碳交易制度还可以促进企业之间的碳排放交易，推动碳市场的形成和发展。企业可以根据自身的碳排放情况进行碳排放交易，实现碳排放的灵活配置，降低整体减排成本。这有助于形成一个更加健康和活跃的碳市场，推动碳排放减排工作的深入开展。总的来说，阶梯碳交易制度是一种创新的碳排放交易制度，可以有效地激励企业减少碳排放，促进绿色发展，推动碳市场的发展。通过不断完善和推广这种制度，我们有望实现更加可持续的发展目标，建设更加清洁和美丽的环境。

2 含阶梯式碳交易机制与电制氢的IES 运行框架

集合多种能源形式的IES通过多种能源与供能设备满足内部的能源需求。本文在传统模型的基础上,引入阶梯式碳交易机制,同时细化考虑了 P2G装置两阶段运行过程中氢能的高效利用以及 CHP 设备热电比可调特性,具体框架如图1所示。

3 阶梯式碳交易机制模型

碳交易机制是通过建立合法的碳排放权,并允许生产商到市场进行碳排放权交易进而达到控制碳排放的目的。监管部门首先为每个碳排放源分配碳排放权配额,生产商结合自身的配额进行合理生产与排放。若实际碳排放低于所分配的配额,则可将剩余配额参与到碳交易市场进行交易,反之则需要购买碳排放权配额。阶梯式碳交易机制模型主要包含碳排放权配额模型、实际碳排放模型、阶梯式碳排放交易模型。

阶梯式定价机制划分了多个购买区间,随着需要购买的碳排放权配额越多,相应区间的购价越高。阶梯式碳交易成本为:

4 阶梯式碳交易机制模型

为验证所提低碳经济调度策略的有效性,设置算例进行验证。以一天 24 h 为一个周期进行优化调度。

5 程序运行结果

1）电功率平衡

2）热功率平衡

3）氢气平衡

3）天然气平衡

6 matlab程序

clc
clear
close all

%% 决策变量初始化
P_CHP_e=sdpvar(1,24); %CHP的输出电功率
P_CHP_h=sdpvar(1,24); %CHP的输出热功率
P_g_CHP=sdpvar(1,24); %CHP消耗天然气功率
P_e_EL=sdpvar(1,24);  %EL设备的耗电量
P_EL_H=sdpvar(1,24);  %EL电解槽的产氢功率
P_H_MR=sdpvar(1,24);  %输入MR设备的氢能功率
P_MR_g=sdpvar(1,24);  %MR设备输出的天然气功率
P_H_HFC=sdpvar(1,24); %输入HFC设备的氢能功率
P_HFC_e=sdpvar(1,24); %HFC设备输出的电功率 
P_HFC_h=sdpvar(1,24); %HFC设备输出的热功率
P_DG=sdpvar(1,24); %风电消纳功率
P_g_GB=sdpvar(1,24); %输入GB设备的天然气功率
P_GB_h=sdpvar(1,24); %GB设备输出的热功率
%储能部分(电ES1、热ES2、气ES3、氢ES4)
P_ES1_cha=sdpvar(1,24);P_ES2_cha=sdpvar(1,24);P_ES3_cha=sdpvar(1,24);P_ES4_cha=sdpvar(1,24); %充放功率
P_ES1_dis=sdpvar(1,24);P_ES2_dis=sdpvar(1,24);P_ES3_dis=sdpvar(1,24);P_ES4_dis=sdpvar(1,24);
S_1=sdpvar(1,24);S_2=sdpvar(1,24);S_3=sdpvar(1,24);S_4=sdpvar(1,24); %各储能的实时容量状态
%引入充放标志二进制变量
B_ES1_cha=binvar(1,24);B_ES2_cha=binvar(1,24);B_ES3_cha=binvar(1,24);B_ES4_cha=binvar(1,24); %充标志
B_ES1_dis=binvar(1,24);B_ES2_dis=binvar(1,24);B_ES3_dis=binvar(1,24);B_ES4_dis=binvar(1,24); %放标志
P_e_buy=sdpvar(1,24); %购电功率
P_g_buy=sdpvar(1,24); %购气功率
ww=275*rand(1,24);
%% 导入常数参数
%电负荷
P_e_load=1.1*[717.451523545706;695.290858725762;689.750692520776;698.060941828255;745.152354570637;808.864265927978;836.565096952909;872.576177285319;886.426592797784;900.277008310249;894.736842105263;883.656509695291;875.346260387812;864.265927977839;864.265927977839;868.421052631579;876.731301939058;889.196675900277;880.886426592798;864.265927977839;836.565096952909;817.174515235457;772.853185595568;745.152354570637]';
%热负荷
P_h_load=0.9*[864.265927977839;941.828254847645;958.448753462604;955.678670360111;988.919667590028;997.229916897507;903.047091412742;833.795013850416;786.703601108033;703.601108033241;664.819944598338;626.038781163435;595.567867036011;590.027700831025;565.096952908587;639.889196675900;714.681440443213;806.094182825485;811.634349030471;831.024930747922;811.634349030471;808.864265927978;800.554016620499;808.864265927978]';
%天然气负荷
P_g_load=[229.916897506925;224.376731301939;216.066481994460;221.606648199446;224.376731301939;252.077562326870;268.698060941828;288.088642659280;299.168975069252;288.088642659280;293.628808864266;282.548476454294;279.778393351801;271.468144044321;271.468144044321;268.698060941828;277.008310249307;293.628808864266;307.479224376731;304.709141274238;293.628808864266;285.318559556787;277.008310249307;265.927977839335]';
%风电预测出力
P_DG_max=1.2*[850.415512465374;864.265927977839;886.426592797784;891.966759002770;894.736842105263;849.030470914127;833.795013850416;653.739612188366;556.786703601108;501.385041551247;432.132963988920;310.249307479224;240.997229916897;252.077562326870;265.927977839335;296.398891966759;343.490304709141;354.570637119114;426.592797783934;526.315789473684;675.900277008310;742.382271468144;854.570637119114;878.116343490305]'+ww;
%购电价格
c_e_buy=1.2*[0.38*ones(1,7),0.68*ones(1,4),1.2*ones(1,3),0.68*ones(1,4),1.2*ones(1,4),0.38*ones(1,2)];
%购气价格
c_g_buy=0.45*ones(1,24);
%% 导入约束条件
C=[];
C=[C,
   P_CHP_e==0.92*P_g_CHP, %CHP的电-气能量转换约束
   0<=P_g_CHP<=600, %CHP消耗的气功率上下限约束
   0.5*P_CHP_e<=P_CHP_h, %热电比上下限
   P_CHP_h<=2.1*P_CHP_e, %热电比上下限
   -0.2*600<=P_g_CHP(2:24)-P_g_CHP(1:23)<=0.2*600, %CHP的爬坡约束(1-24时段)
  ];
C=[C,
   P_EL_H==0.88*P_e_EL, %EL(电解槽)的氢-电能量转换约束
   0<=P_e_EL<=500, %EL的消耗电功率的上下限约束
   -0.2*500<=P_e_EL(2:24)-P_e_EL(1:23)<=0.2*500, %EL的爬坡约束(1-24时段)
  ];
C=[C,
   P_MR_g==0.6*P_H_MR, %MR(甲烷反应器)的气-氢能量转换约束
   0<=P_H_MR<=250, %MR消耗的氢功率的上下限约束
   -0.2*250<=P_H_MR(2:24)-P_H_MR(1:23)<=0.2*250, %MR的爬坡约束(1-24时段)
  ];
C=[C,
   P_HFC_e==0.85*P_H_HFC, %HFC(氢燃料电池)的电-氢能量转换约束
   0<=P_H_HFC<=250, %HFC消耗的氢功率上下限约束
   0.5*P_HFC_e<=P_HFC_h, %HFC的热电比上下限
   P_HFC_h<=2.1*P_HFC_e, %HFC的热电比上下限
   -0.2*250<=P_H_HFC(2:24)-P_H_HFC(1:23)<=0.2*250, %HFC的爬坡约束(1-24时段)
  ];
C=[C,
   0<=P_DG<=P_DG_max, %风电出力约束
   P_GB_h==0.95*P_g_GB, %GB的热-气能量转换约束
   0<=P_g_GB<=800, %GB的出力上下限约束
   -0.2*800<=P_g_GB(2:24)-P_g_GB(1:23)<=0.2*800, %GB的爬坡约束(1-24时段)
  ];
C=[C,
   0<=P_ES1_cha<=B_ES1_cha*0.5*450, %储电设备的最大充电功率约束
   0<=P_ES2_cha<=B_ES2_cha*0.5*500, %储热设备的最大充热功率约束
   0<=P_ES3_cha<=B_ES3_cha*0.5*150, %储气设备的最大充气功率约束
   0<=P_ES4_cha<=B_ES4_cha*0.5*200, %储氢设备的最大充氢功率约束
   0<=P_ES1_dis<=B_ES1_dis*0.5*450, %储电设备的最大放电功率约束
   0<=P_ES2_dis<=B_ES2_dis*0.5*500, %储热设备的最大放热功率约束
   0<=P_ES3_dis<=B_ES3_dis*0.5*150, %储气设备的最大放气功率约束
   0<=P_ES4_dis<=B_ES4_dis*0.5*200, %储氢设备的最大放氢功率约束
   S_1(1)==0.3*450, %储电设备的初始容量
   S_2(1)==0.3*500, %储热设备的初始容量
   S_3(1)==0.3*150, %储气设备的初始容量
   S_4(1)==0.3*200, %储氢设备的初始容量
   %始末状态守恒约束
   S_1(24)==S_1(1), 
   S_2(24)==S_2(1),
   S_3(24)==S_3(1),
   S_4(24)==S_4(1),   
   %充放状态唯一
   B_ES1_cha+B_ES1_dis<=1, 
   B_ES2_cha+B_ES2_dis<=1,
   B_ES3_cha+B_ES3_dis<=1,
   B_ES4_cha+B_ES4_dis<=1,
   %储能容量上下限约束
   0.2*450<=S_1<=0.9*450,
   0.2*500<=S_2<=0.9*500,
   0.2*150<=S_3<=0.9*150,
   0.2*200<=S_4<=0.9*200,
   %储能容量变化约束
   S_1(2:24)==S_1(1:23)+0.95*P_ES1_cha(2:24)-P_ES1_dis(2:24)/0.95,
   S_2(2:24)==S_2(1:23)+0.95*P_ES2_cha(2:24)-P_ES2_dis(2:24)/0.95,
   S_3(2:24)==S_3(1:23)+0.95*P_ES3_cha(2:24)-P_ES3_dis(2:24)/0.95,
   S_4(2:24)==S_4(1:23)+0.95*P_ES4_cha(2:24)-P_ES4_dis(2:24)/0.95,
  ];
C=[C,
   P_e_buy==P_e_load+P_e_EL+P_ES1_cha-P_ES1_dis-P_DG-P_CHP_e-P_HFC_e, %电功率平衡约束
   P_HFC_h+P_CHP_h+P_GB_h==P_h_load+P_ES2_cha-P_ES2_dis, %热功率平衡约束
   P_g_buy==P_g_load+P_ES3_cha-P_ES3_dis+P_g_CHP+P_g_GB-P_MR_g, %气功率平衡约束
   P_EL_H==P_H_MR+P_H_HFC+P_ES4_cha-P_ES4_dis, %氢功率平衡约束
   0<=P_e_buy<=5000, %购电功率约束
   0<=P_g_buy<=5000, %购气功率约束
  ];
%% 导入目标函数
%公式5 碳排放权配额模型
E_e_buy=0.68*sum(P_e_buy); %购电配额
E_CHP=0.102*3.6*sum(P_CHP_h+6/3.6*P_CHP_e); %CHP配额
E_GB=0.102*3.6*sum(P_GB_h); %GB配额
E_IES=E_e_buy+E_CHP+E_GB; %IES总碳排放配额
%公式6 实际碳排放模型
E_e_buy_a=1.08*sum(P_e_buy);
E_CHP_a=0.065*3.6*sum(P_CHP_h+6/3.6*P_CHP_e); 
E_GB_a=0.065*3.6*sum(P_GB_h); 
E_MR_a=1*sum(P_MR_g); %实际MR减少的碳排放
E_IES_a=E_e_buy_a+E_CHP_a+E_GB_a;
E=E_IES_a-E_MR_a; %实际IES总碳排放
%阶梯碳交易成本(分段线性化)
E_v=sdpvar(1,5)%每段区间内的长度,分为5段,每段长度是2000
lamda=0.250;%碳交易基价
C=[C,
   E==sum(E_v),%总长度等于E
   0<=E_v(1:4)<=2000,%除了最后一段，每段区间长度小于等于2000
   0<=E_v(5),
  ];
%碳交易成本
C_CO2=0;
for v=1:5
    C_CO2=C_CO2+(lamda+(v-1)*0.25*lamda)*E_v(v);
end
%总目标函数
Cost=0;
for t=1:24
    Cost=Cost+c_e_buy(t)*P_e_buy(t)'+c_g_buy(t)*P_g_buy(t)'+C_CO2+0.13*sum(P_DG_max(t)-P_DG(t)); %购能成本+碳交易成本+弃风成本
end
%% 求解器相关配置
ops=sdpsettings('solver','cplex','verbose',2,'usex0',0);
ops.cplex.mip.tolerances.mipgap=1e-6;
%% 进行求解计算         
result=optimize(C,Cost,ops);
if result.problem==0 
else
    error('求解出错');
end  
Cost=value( Cost);
display(['通过Yalmip求得的最优规划值为 : ', num2str(Cost)]);

%% 画图输出运行结果

%电功率情况
Plot_e=[];
for t=1:24
    Plot_e(1,t)=P_e_buy(t);
    Plot_e(2,t)=P_HFC_e(t);
    Plot_e(3,t)=P_ES1_dis(t);
    Plot_e(4,t)=P_DG(t);
    Plot_e(5,t)=P_CHP_e(t);
    Plot_e(6,t)=P_e_EL(t);
    Plot_e(7,t)=P_ES1_cha(t);
end

figure
bar(Plot_e(1:5,:)','stacked');
hold on
bar(-Plot_e(6:7,:)','stacked');
hold on
plot(P_e_load,'r-o','LineWidth',1.5);
xlabel('时间/h');
ylabel('功率/kW');
title('电功率平衡图');
legend('购电','HFC','放电','风电','CHP','EL','储电','电负荷');

%热功率情况
Plot_h=[];
for t=1:24
    Plot_h(1,t)=P_HFC_h(t);
    Plot_h(2,t)=P_CHP_h(t);
    Plot_h(3,t)=P_GB_h(t);
    Plot_h(4,t)=P_ES2_dis(t);
    Plot_h(5,t)=P_ES2_cha(t);
end

figure
bar(Plot_h(1:4,:)','stacked');
hold on
bar(-Plot_h(5,:)','stacked');
hold on
plot(P_h_load,'k-o','LineWidth',1.5);
xlabel('时间/h');
ylabel('功率/kW');
title('热功率平衡图');
hold on
legend('HFC','CHP','GB','热放','热储','热负荷');

%氢功率情况
Plot_H=[];
for t=1:24
    Plot_H(1,t)=P_EL_H(t);
    Plot_H(2,t)=P_ES4_dis(t);
    Plot_H(3,t)=P_H_MR(t);
    Plot_H(4,t)=P_H_HFC(t);
    Plot_H(5,t)=P_ES4_cha(t);
end

figure
bar(Plot_H(1:2,:)','stacked');
hold on
bar(-Plot_H(3:5,:)','stacked');
xlabel('时间/h');
ylabel('功率/kW');
title('氢功率平衡图');
legend('P2G','氢放','MR','HFC','氢储');


%气功率情况
Plot_g=[];
for t=1:24
    Plot_g(1,t)=P_g_buy(t);
    Plot_g(2,t)=P_ES3_dis(t);
    Plot_g(3,t)=P_MR_g(t);
    Plot_g(4,t)=P_g_GB(t);
    Plot_g(5,t)=P_g_CHP(t);
    Plot_g(6,t)=P_ES3_cha(t);
end

figure
bar(Plot_g(1:3,:)','stacked');
hold on
bar(-Plot_g(4:6,:)','stacked');
hold on
plot(P_g_load,'k-*','LineWidth',1.5);
xlabel('时间/h');
ylabel('功率/kW');
title('气负荷功率平衡图');
legend('购气','放气','甲烷反应器','GB','CHP','储气','气负荷');
P_DG=value(P_DG);