基于共享储能电站的工业用户日前优化经济调度MATLAB（yalmip+cplex）

1. 文章介绍

摘要：文章提出一种基于共享储能电站的工业用户日前优化经济调度方法。首先提出共享储能电站的概念,分析其商业运营模式。然后将共享储能电站应用到工业用户经济优化调度中,通过协调各用户使用共享储能电站进行充电和放电的功率,实现用户群日运行成本最优。最后以江苏省 3 个工业用户进行算例仿真,与用户不配置储能和用户独立配置储能场景对比,得出引入共享储能电站可以显著降低用户群日运行成本,并对储能电站年服务费收益、静态投资回收年限和投资回报率与共享储能电站服务费定价间的关系做进一步的研究。

1.1 共享储能电站的运营模式

        共享储能服务即用户使用共享储能电站的储能装置满足用户内部充电和放电的需求,不受时间和容量的限制。用户通过向储能电站运营商缴纳服务费用换取共享储能服务。储能服务费用是用户使用共享储能电站充电和放电,需向共享储能电站缴纳的费用。

1.2 基于共享储能的优化调度模型

1.2.1 目标函数

        接入共享储能电站的用户群以典型日运行成本最优为目标函数,确定储能电站的容量、最大充放电功率和各个时段的充放电功率、用户各个时段利用储能电站进行充电和放电的功率。用户群的日运行成本包括从电网购电费用和向储能电站缴纳的服务费。

1.2.2 约束条件

• 电功率平衡约束

• 用户使用共享储能电站充放电功率约束

• 储能电站荷电状态连续性约束

• 储能电站充放电功率约束

• 储能电站充放电功率平衡约束

1.3 求解方法:大M法

        约束条件储能电站充放电功率约束为非线性约束,  采用 Big-M 法对非线性约束进行线性化,转化为约束条件式( 9) ,其中 M 为足够大的常数。

2.  算例分析，场景对比

场景 1: 各用户不配置储能;  场景 2: 各用户内部独立配置储能;  场景 3: 各用户接入共享储能电站。

3. 代码

%% 复现场景3
clc
clear
close all

%% 参数设置
P_grid=sdpvar(3,24); %用户从电网购买的电功率
P_ess_b=sdpvar(3,24); %用户使用共享储能电站放电的功率
P_ess_s=sdpvar(3,24); %用户使用共享储能电站充电的功率
U_ess_b=binvar(3,24); %用户使用共享储能的放电状态位,取1时为放电,0为未放电
U_ess_s=binvar(3,24); %用户使用共享储能的充电状态位,取1时为充电,0为未充电
E=sdpvar(1,24); %储能电站的荷电状态
P_abs=sdpvar(1,24); %储能电站的充电功率
P_relea=sdpvar(1,24); %储能电站的放电功率
E_init=sdpvar(1,1); %储能电站的初始容量
E_max=sdpvar(1,1); %储能电站的预配置容量
P_max=sdpvar(1,1); %储能电站的预配置最大充放电功率
U_abs=binvar(1,24); %储能电站的放电状态位,取1时为放电,0为未放电
U_relea=binvar(1,24); %储能电站的充电状态位,取1时为充电,0为未充电

%% 用户参数导入
gamma=[0.37*ones(1,8),1.36*ones(1,4),0.82*ones(1,5),1.36*ones(1,4),0.82*ones(1,3)]; %电网的"峰-平-谷"电价
delta=0.33*ones(1,24); %共享储能电站的服务费
P_load(1,:)=[80,75,75,75,75,75,75,75,95,155,160,180,160,165,150,165,160,175,160,130,140,100,90,85];
P_load(2,:)=[200,50,150,190,25,50,30,20,40,170,220,150,65,160,170,150,230,125,140,200,125,70,110,65];
P_load(3,:)=[105,100,95,100,100,105,110,110,120,110,115,120,135,110,120,125,130,125,115,130,120,120,115,120];
P_pv(1,:)=[0,0,0,0,0,0,0,0,75,275,540,675,750,750,330,170,75,10,0,0,0,0,0,0];
P_pv(2,:)=zeros(1,24);P_pv(3,:)=zeros(1,24);
P_wind(1,:)=zeros(1,24);
P_wind(2,:)=[150,150,175,165,140,130,120,90,50,55,80,100,135,135,130,110,60,50,50,70,90,120,140,120];
P_wind(3,:)=[130,145,155,135,120,125,115,95,40,30,50,85,110,115,110,70,35,40,40,50,80,80,115,95];

%% 约束条件
C=[]; %约束条件矩阵初始化
for i=1:3
    for t=1:24
        C=[C,
           P_pv(i,t)+P_wind(i,t)+P_grid(i,t)+P_ess_b(i,t)-P_ess_s(i,t)-P_load(i,t)==0, %电功率平衡约束
           0<=P_ess_b(i,t)<=1000*U_ess_b(i,t), %用户使用共享储能电站放电功率约束
           0<=P_ess_s(i,t)<=1000*U_ess_s(i,t), %用户使用共享储能电站充电功率约束
           U_ess_b(i,t)+U_ess_s(i,t)<=1, %避免同时出现充放电的情况,故增加此约束
           P_grid(i,t)>=0,
          ]; %电功率平衡约束
    end
end
%储能电站荷电状态连续性约束
C=[C,E(1)==E_init+0.95*P_abs(1)-P_relea(1)/0.95,]; %1时段和其他时段不同,涉及E(0)
for t=2:24
    C=[C,E(t)==E(t-1)+0.95*P_abs(t)-P_relea(t)/0.95,]; 
end
for t=1:24
    C=[C,0.1*E_max<=E(t)<=0.9*E_max,];
end
C=[C,E_init==0.2*E_max,E(24)==E_init,];
%储能电站的充放电功率约束,用文章2.3中的Big-M法进行线性化处理
M=1E8; %这里的M是个很大的数
for t=1:24
    C=[C,
       0<=P_abs(t)<=P_max,
       0<=P_abs(t)<=U_abs(t)*M,
       0<=P_relea(t)<=P_max,
       0<=P_relea(t)<=U_relea(t)*M,
       U_abs(t)+U_relea(t)<=1,
      ];
end
%储能电站充放电功率平衡约束
for t=1:24
    C=[C,
       (P_ess_b(1,t)-P_ess_s(1,t))+(P_ess_b(2,t)-P_ess_s(2,t))+(P_ess_b(3,t)-P_ess_s(3,t))==P_relea(t)-P_abs(t),
      ];
end

%% 目标函数
C1=sum(gamma.*P_grid(1,:))+sum(gamma.*P_grid(2,:))+sum(gamma.*P_grid(3,:)); %用户群从电网购电费用
C2=sum(delta.*(P_ess_b(1,:)+P_ess_s(1,:)))+sum(delta.*(P_ess_b(2,:)+P_ess_s(2,:)))+sum(delta.*(P_ess_b(3,:)+P_ess_s(3,:))); %用户群向储能电站缴纳的服务费
F=C1+C2; %总目标函数
%% 模型求解
ops=sdpsettings('solver','cplex','verbose',2,'usex0',0);
ops.cplex.mip.tolerances.mipgap=1e-6;
result=optimize(C,F,ops);

if result.problem == 0 % problem =0 代表求解成功
else
    error('求解出错');
end

%% 数据导出及画图
E_max=value(E_max);P_max=value(P_max);
display(['最优储能电站容量规划值为 : ', num2str(E_max),' kWh']);
display(['最优储能电站充放电功率最大值为 : ', num2str(P_max),' kW']);
P_grid=value(P_grid);P_ess_b=value(P_ess_b);P_ess_s=value(P_ess_s);
E=value(E);P_abs=value(P_abs);P_relea=value(P_relea);
%用户A的电负荷平衡曲线
for t=1:24
    if P_grid(1,t)<0.01
        P_grid(1,t)=0; %由于Cplex求解精度,把极小的值设为0
    end
end
for t=1:24
    Plot_UserA(1,t)=P_grid(1,t);
    Plot_UserA(2,t)=P_ess_b(1,t);
    Plot_UserA(3,t)=-1*P_ess_s(1,t);
end
Plot_UserA=Plot_UserA';
figure
bar(Plot_UserA,'stacked');
hold on
plot(P_load(1,:),'k-o','LineWidth',1.5);
hold on
plot(P_pv(1,:),'b-o','LineWidth',1.5);
xlabel('时间/h');
ylabel('功率/kW');
title('用户A的电负荷平衡曲线');
legend('从电网购电功率','使用储能电站的放电功率','使用储能电站的充电功率','用户A电负荷','光伏出力');
legend('boxoff');
box off
%用户B的电负荷平衡曲线
for t=1:24
    if P_grid(2,t)<0.01
        P_grid(2,t)=0; %由于Cplex求解精度,把极小的值设为0
    end
end
for t=1:24
    Plot_UserB(1,t)=P_grid(2,t);
    Plot_UserB(2,t)=P_ess_b(2,t);
    Plot_UserB(3,t)=-1*P_ess_s(2,t);
end
Plot_UserB=Plot_UserB';
figure
bar(Plot_UserB,'stacked');
hold on
plot(P_load(2,:),'k-o','LineWidth',1.5);
hold on
plot(P_wind(2,:),'b-o','LineWidth',1.5);
xlabel('时间/h');
ylabel('功率/kW');
title('用户B的电负荷平衡曲线');
legend('从电网购电功率','使用储能电站的放电功率','使用储能电站的充电功率','用户B电负荷','风电出力');
legend('boxoff');
box off
%用户C的电负荷平衡曲线
for t=1:24
    if P_grid(3,t)<0.01
        P_grid(3,t)=0; %由于Cplex求解精度,把极小的值设为0
    end
end
for t=1:24
    Plot_UserC(1,t)=P_grid(3,t);
    Plot_UserC(2,t)=P_ess_b(3,t);
    Plot_UserC(3,t)=-1*P_ess_s(3,t);
end
Plot_UserC=Plot_UserC';
figure
bar(Plot_UserC,'stacked');
hold on
plot(P_load(3,:),'k-o','LineWidth',1.5);
hold on
plot(P_wind(3,:),'b-o','LineWidth',1.5);
xlabel('时间/h');
ylabel('功率/kW');
title('用户C的电负荷平衡曲线');
legend('从电网购电功率','使用储能电站的放电功率','使用储能电站的充电功率','用户C电负荷','风电出力');
legend('boxoff');
box off
%共享储能电站充放电功率和电量状态曲线
for t=1:24
    Plot_Station(1,t)=-1*P_abs(t);
    Plot_Station(2,t)=P_relea(t);
end
Plot_Station=Plot_Station';
EE=E';
figure
yyaxis left; % 激活左边的轴
bar(Plot_Station,'stacked');
title('共享储能电站充放电功率和电量状态曲线');
xlabel('时间/h');
ylabel('功率/kW'); % 给左y轴添加轴标签
yyaxis right; % 激活右边的轴
plot(EE,'r-o','LineWidth',1.5);
ylabel('电量/kWh'); % 给右y轴添加轴标签
legend('储能电站充电功率','储能电站放电功率','储能电站电量');
legend('boxoff');
grid on
box off     

4. 程序运行结果