基于多能互补的热电联供型微网优化运行（matlab代码）

目录

1 主要内容

多能互补模型

算例分析

2 部分代码

3 程序结果

4 下载链接

---

1 主要内容

该程序基本复现《基于多能互补的热电联供型微网优化运行》，在需求侧对负荷类型进行分类，利用电负荷的弹性和系统供热方式的多样性，构建含电负荷时移、削减响应及热负荷供能方式响应的综合能源需求响应模型，并提出响应补偿机制。在此基础上，以系统运行成本与响应补偿成本之和最小为目标，综合考虑供需双侧设备运行和可调度负荷资源约束，建立基于多能互补的 CHP-MG 优化运行数学模型。同时，为了体现模型的准确性，代码对比了热负荷参与、电负荷参与以及电热负荷均参与或者均不参与四种常见下的模型调度结果，体现了所构模型的经济型。

• 多能互补模型

• 算例分析

2 部分代码

%%  赋值
MT=intvar(1,24,'full');
​
Pgrid=intvar(1,24,'full');  %  购 、 售
​
Pnas=intvar(2,24,'full'); % 充、 放
​
​
H=intvar(1,24,'full');%锅炉
Q=intvar(1,24,'full');%可中断电负荷
​
Hti=intvar(1,24,'full');%充热
Hto=intvar(1,24,'full');%放热
xx=intvar(1,24,'full');%%%%%时平负荷量
yy=intvar(1,24,'full');
​
%% 0-1赋值
I_Hti=binvar(1,24,'full');%充热
I_Hto=binvar(1,24,'full');%放热
I_MT=binvar(1,24,'full');
I_Pnas=binvar(2,24,'full');% 1运行 0停止
I_Q=binvar(1,24,'full');
%%   目标函数
for i=1:24%发电成本
    Cf(1,i)=Cgas(i)*(aF(1)*MT(1,i)+bF(1)*I_MT(1,i));
end
for i=1:24%%%%余热回收
    H_cycle(1,i)=aH(1)*MT(1,i)+bH(1)*I_MT(1,i);
end
for k=1:24 %PCC交互成本 % 1-5，23-24 谷 % 6-12，19-22 峰 % 13-18 平 
    if k>=1&&k<7
        Cgrid(1,k)=Pgrid(1,k).*buy(3);
    elseif k>=7&&k<13
        Cgrid(1,k)=Pgrid(1,k).*buy(1);
    elseif k>=13&&k<19
        Cgrid(1,k)=Pgrid(1,k).*buy(2);
    elseif k>=19&&k<23
        Cgrid(1,k)=Pgrid(1,k).*buy(1);
    else
        Cgrid(1,k)=Pgrid(1,k).*buy(3);
    end
end
for k=1:24 %需求响应单位成本
   if k>=1&&k<7
        bu_q(1,k)=0.9*buy(1);
        bu_x(1,k)=0.5*buy(1);
        bu_p(1,k)=0.6*Cgas(1,k);
    elseif k>=7&&k<13
        bu_q(1,k)=0.9*buy(1);
        bu_x(1,k)=0.5*buy(1);
        bu_p(1,k)=0.6*Cgas(1,k);
    elseif k>=13&&k<19
        bu_q(1,k)=0.9*buy(1);
        bu_x(1,k)=0.5*buy(1);
        bu_p(1,k)=0.6*Cgas(1,k);
    elseif k>=19&&k<23
        bu_q(1,k)=0.9*buy(1);
        bu_x(1,k)=0.5*buy(1);
        bu_p(1,k)=0.6*Cgas(1,k);
    else
        bu_q(1,k)=0.8*buy(1);
        bu_x(1,k)=buy(1)/2;
        bu_p(1,k)=0.6*Cgas(k);
    end
end
for k=1:24  %% 切除负荷成本
        if k>=7&&k<=12
          Ck1(1,k)=(Q(1,k).*bu_q(1,k))+xx(1,k)*bu_x(1,k);
    elseif  k>=19&&k<=20
          Ck1(1,k)=(Q(1,k).*bu_q(1,k))+xx(1,k)*bu_x(1,k);
       else
          Ck1(1,k)=(Q(1,k).*bu_q(1,k))+xx(1,k)*bu_x(1,k);
        end
end
for i=1:24%锅炉成本
    Ch(1,i)=Cgas(i)*(H(1,i))/LHV;
end
F=0;%目标函数
mm=3.1;
for k=1:24 %1.8
     F=F+Cf(1,k)+Cgrid(1,k)+Ch(1,k)+(Pnas(1,k)+Pnas(2,k))*0.024;
end
for k=1:24 %SOC值
    SOC(k)=(500+sum(Pnas(1,1:k).*I_Pnas(1,1:k)-(Pnas(2,1:k)).*I_Pnas(2,1:k)))/1000;   
end
begin=500;
% for i=1:24%%热储能
%     L(1,i)=begin*h_n+h_charge*Hti(1,i)-Hto(1,i);%%%热储能容量
%     begin=L(1,i);
% end
​
%%  约束条件
constraints=[];
%% 状态约束
for k=1:24  %Pgrid状态 %Pnas状态
    constraints=[constraints,I_Pnas(1,k)+I_Pnas(2,k)<=1];
%     constraints=[constraints,I_Hti(1,k)+I_Hto(1,k)<=1];
end
constraints=[constraints,sum(I_Pnas(1,1:24)+I_Pnas(2,1:24))<=14];
%% 上下限约束
for k=1:24     
    constraints=[constraints,25.*I_MT(1,k)<=MT(1,k)<=145.*I_MT(1,k)]; 
    
    constraints=[constraints,Pgrid_min<=Pgrid(1,k)<=Pgrid_max]; 
​
    constraints=[constraints,Pnas_min.*I_Pnas(1,k)<=Pnas(1,k)<=Pnas_max.*I_Pnas(1,k)]; 
    constraints=[constraints,Pnas_min.*I_Pnas(2,k)<=Pnas(2,k)<=Pnas_max.*I_Pnas(2,k)];    
end
%%% MT爬坡率
for i=1:23
    constraints=[constraints,-55<=MT(1,i+1)-MT(1,i)<=65];
end
%   PCC深度限制
 for k=1:23
     constraints=[constraints,-90<=Pgrid(1,k+1)-Pgrid(1,k)<=90];
 end
 %%荷电状态
 for k=1:24
      constraints=[constraints,SOC_min<=300+sum(Pnas(1,1:k)-Pnas(2,1:k))<=SOC_max];      
 end
 %%储能充放电深度限制
 for k=1:23
     constraints=[constraints,-50<=Pnas(1,k+1)-Pnas(2,k+1)-Pnas(1,k)+Pnas(2,k)<=50];
 end
     constraints=[constraints,sum(Pnas(1,1:24))==sum(Pnas(2,1:24))];
%%%锅炉上下限、爬坡率
for i=1:24
    constraints=[constraints,30<=H(1,i)<=H_max];
end
for i=1:23   
   constraints=[constraints,-90<=H(1,i+1)-H(1,i)<=90];
end
​
%%功率平衡
for k=1:24
constraints=[constraints,MT(1,k)-Pnas(1,k)+Pnas(2,k)+Pgrid(1,k)+Pwind(1,k)==load(k)+L5(1,k)+L6(1,k)+L7(1,k)];
end

3 程序结果