基于EFC、ELCC定义下的风电可信容量计算matlab程序

%% 4.1计算风电等值容量（EFC定义下）

%相关可靠性计算方法参见百度文库“含风电的发电系统可靠性评估（MC法）matlab程序”clear

load process

loadresultwind

load result0

% 计算EFC定义下风电等值容量

P0=0;%设定初值Pe下限为0

R0=EENS2;

P=[];%设定变量P、R存储每次迭代后Pe和对应的可靠性指标

R=[];

% 计算初值P(1)的可靠性指标,P(1)对应等值容量上限，和风电容量相同

P(1)=Generator.Wind.Pr*N;

Generator.Norm=[Generator.Norm;[P(1),0]];%第12行【15 0】是EFC定义下的理想机组

I=0;

sumDNS=0;

DNS=zeros(100000,1);

K=rand(100000,12);%1-11常规12理想机组

for k1=1:100000

Pout=zeros(12,1);

%得到一次抽样常规机组状态

fori=1:12

if K(k1,i)>Generator.Norm(i,2)

Pout(i)=Generator.Norm(i,1);

else

Pout(i)=0;

end

end

sumPout=sum(Pout);%发电系统出力

% 统计系统状态参数

% 持续峰荷185MW

RBTSload=185;

DNS(k1)=max([0;RBTSload-sumPout]);

if DNS(k1)>0

I=I+1; %I用来记录发生却负荷的次数

sumDNS=sumDNS+DNS(k1);

end

end

% 计算指标

EENS=sumDNS/100000*8760;

R(1)=EENS;

% 以下进行可信容量的搜索（弦截法）

g=1;

deta=10;

while abs(EENS1-R(g))>=deta

P(g+1)=(EENS1-R0)*P(g)/(R(g)-R0);

Generator.Norm=[Generator.Norm(1:11,:);[P(g+1),0]];

% 重复进行含理想机组的可靠性评估

I=0;

sumDNS=0;

DNS=zeros(100000,1);

K=rand(100000,12);%1-11常规12理想机组

for k1=1:100000

Pout=zeros(12,1);

%得到一次抽样常规机组状态

fori=1:12

if K(k1,i)>Generator.Norm(i,2)

Pout(i)=Generator.Norm(i,1);

else

Pout(i)=0;

end

end

sumPout=sum(Pout);%发电系统出力

% 统计系统状态参数

% 持续峰荷185MW

RBTSload=185;

DNS(k1)=max([0;RBTSload-sumPout]);

if DNS(k1)>0

I=I+1; %I用来记录发生却负荷的次数sumDNS=sumDNS+DNS(k1);

end

end

EENS=sumDNS/100000*8760;

R(g+1)=EENS;

g=g+1;

end

Peq=P(g);

rate_credict=Peq/15;

disp('风电等值容量');

disp(Peq);

disp('风电容量可信度');

disp(rate_credict);

resultEFC.Peq=Peq;

resultEFC.rate_credict=rate_credict;

save('result_WindCredict.mat','resultEFC') %% 4.2计算风电可信容量（ELCC定义下）clear

load result0

load process

loadresultwind

D=[];

R=[];

% 计算200MW负荷对应的可靠性指标R1 I=0;

sumDNS=0;

DNS=zeros(200000,1);

K=rand(200000,12);%1-11常规12风电pwind=zeros(200000,1);

for k=1:200000

Pout=zeros(12,1);

%得到一次抽样常规机组状态

fori=1:11

if K(k,i)>Generator.Norm(i,2)

Pout(i)=Generator.Norm(i,1);

else

Pout(i)=0;

end

end

%得到一次抽样风电场的出力

t=0;

fori=1:M

t=t+StateFORWeibull6(i,2);

if K(k,12)-t>=0

continue

else

Pout(12)=StateFORWeibull6(i,1);

pwind(k)= Pout(12);

break