togrid.m

基于长期负荷预测的线性规划燃机出力能量优化控制

clear;
format long;
%STA=[0 864 1512 2352 3024 3864 4536 5376 6048 6888 7560 8400];
sta=6048;%(31+4+28+31+4+28+35+28+35+28+35+28+35)*24;

[E_T,E_D,E_O,Q_T,Q_D,Q_O]=creat(sta);
Le=24;%length(E_T);
E_gt=zeros(Le,1);
P_ep=zeros(Le,1);
Q_ac=zeros(Le,1);
PER=zeros(Le,1);
P_air1=zeros(Le,1);
P_air2=zeros(Le,1);
Obj_Val=zeros(Le,1);
P_gt=zeros(Le,1);
Q_con=zeros(Le,1);
Q1_GT=zeros(Le,1);
Opted=zeros(Le,1);
READED=zeros(Le,1);
FEV1=zeros(Le,1);
FEV2=zeros(Le,1);

t=[1:Le];
p_ng=3.85; p_ele=1.0136;COP_ac=1.4; COP_tes=2.6; COP_dis=5; COP_off=4.82;LHV=41.9; PER_rate=0.75;
a1=0;b1=0.3972;c1=-0.5378;d1=-31.6244;a2=0.000307;b2=0.2418;c2=1.4611;d2=42.784;
d_GT=1;t=35;

a_p=(a1*t+b1);b_p=(c1*t+d1)*d_GT;                 % P_GT=a_p*E_GT+b_p;
a_q=(a2*t+b2);b_q=(c2*t+d2)*d_GT;                 % Q1_GT=a_q*E_GT+b_q;
[E_T,E_D,E_O,Q_T,Q_D,Q_O]=creat(sta);


for i=[1:Le]
    E_tes=E_T(i); E_dis=E_D(i); E_off=E_O(i);% P_gt=0.4*E_gt-40;   %代入可得     0.4*E_gt+P_ep-E_re_dis-E_re_off=P_loads+E_re_tes+40
    Q_tes=Q_T(i); Q_dis=Q_D(i); Q_off=Q_O(i); % Q_gt=0.3*E_gt+40;   %代入可得    E_re_dis*COP_dis+E_re_off*COP_off+Q_ac=Q_dis+Q_off
    [X,EXITF,OPT,FEV1,FEV2]=caculating(a_p,a_q,b_p,b_q,i,p_ng,p_ele,COP_ac,COP_tes,COP_dis,COP_off,LHV,PER_rate,sta,FEV1,FEV2,E_tes,E_dis,E_off,Q_tes,Q_dis,Q_off);
    if EXITF(1)>=0&EXITF(2)>=0
        num=find(OPT==min(OPT));
        if OPT(1)==OPT(2)
            num=1;
        end
        x=X(:,num);
        Obj_Val(i)=min(OPT);
        Opted(i)=EXITF(num);
        READED(i)=1;
    else if EXITF(1)*EXITF(2)<=0
            num1=find(EXITF>=0);
            x=X(:,num1);
            Obj_Val(i)=OPT(num1);
            Opted(i)=EXITF(num1);
            READED(i)=1;
        else if EXITF(1)<0&EXITF(2)<0
                x=zeros(5,1);feval=0;Obj_Val(i)=feval;
                READED(i)=0;
            end
        end
    end
    if READED(i)==1;
        P_gt(i)=a_p*x(2)+b_p;%%%燃机出力
        Q1_GT(i)=a_q*x(2)+b_q;%%%%燃机废热
        E_gt(i)=x(2);%%%%%燃机消耗的燃料热值
        P_ep(i)=x(1);
        P_air1(i)=x(3);%%%%%%新大楼电空调
        P_air2(i)=x(4);
        Q_ac(i)=x(5);%%%%%产生的制冷量
        PER(i)=(Q_ac(i)+P_gt(i))/x(2);
        str=sprintf('E_gt(%2d) = %d kW\t P_ep(%2d) = %d kW\t Obj_Val(%2d) = %d$',i,x(2),i,x(1),i,Obj_Val(i));disp(str);
    else
        P_gt(i)=0;Q1_GT(i)=0; Q_ac(i)=0;PER(i)=0;
        P_air1(i)=Q_D(i)/COP_dis;
        P_ep(i)=E_T(i)+Q_T(i)/COP_tes+E_D(i)*0.5+P_air1(i);
        P_air2(i)=Q_O(i)/COP_off;
    end
    Q_waste(i)=Q1_GT(i)-Q_ac(i)/COP_ac;
    if Q_waste(i)<0
        Q_waste(i)=0;
    end
end

hour=1:24;
figure(1)
plot(hour,P_ep,'-*g');
hold on
plot(hour,P_gt,'-*m');
hold on
P_loads=E_T+E_D/2+Q_T/COP_tes+P_air2;
plot(hour,P_loads,'-*c');
h=legend('P_g_r_i_d','P_M_T','P_l_o_a_d_s',2);
grid on
set(gca,'XTick',1:24)
xlabel('hour');
ylabel('kW');
title('电负荷平衡图')

figure(2)
Q_loads=Q_D+Q_O;
plot(hour,Q_loads,'-*c');
hold on
plot(hour,Q_ac,'-*m')
hold on
Q_con1=P_air1*COP_dis;
Q_con2=P_air2*COP_off;
plot(hour,Q_con1,'-*g')
hold on
plot(hour,Q_con2,'-*r')
hold on
plot(hour,Q_D,'-*b');
grid on
set(gca,'XTick',1:24)
h=legend('Q_l_o_a_d_s','Q_c','Q_c_o_n_1','Q_c_o_n_2','Q_1',2);%Q_con=E_con_D*COP_dis+E_con_O*COP_off表示电空调产生的冷功率大小
xlabel('hour');
ylabel('kW');
title('冷负荷平衡图');

figure(3)
plot(hour,PER,'-*');
hold on
plot(hour,0.75,'-r');
title('PER')
grid on
set(gca,'XTick',1:24)
h=legend('PER_C_C_P',2);
xlabel('hour');
title('效率图');

figure(4)
plot(hour,Obj_Val,'-*b');
hold on
Prefee=(E_T+E_D/2+Q_D/COP_dis+Q_O/COP_off+Q_T/COP_tes)*p_ele;
plot(hour,Prefee,'-*m');
hold on
plot(hour,FEV1,'-*g');
hold on
plot(hour,FEV2,'-*r');
title('Obj_V_a_n')
grid on
set(gca,'XTick',1:24)
h=legend('Obj_V_a_n','Prefee','FEV1','FEV2',2);
xlabel('hour');
ylabel('元');
title('运行费用图');

figure(5)
plot(hour,Q_waste);
grid on
set(gca,'XTick',1:24)
h=legend('Q_w_a_s_t_e_d',2);%Q_con=E_con_D*COP_dis+E_con_O*COP_off表示电空调产生的冷功率大小
xlabel('hour');
ylabel('kW');
title('浪费的燃机废热图');