why my code is not converging... 57 bus load flow

조회 수: 20 (최근 30일)
Ghulam Mohyuddin
Ghulam Mohyuddin 2015년 3월 26일
댓글: SURAJ KUMAR 2023년 6월 11일
clear
basemva = 100; accuracy = 0.001; maxiter = 200;
% IEEE 57-BUS TEST SYSTEM (American Electric Power)
% Bus Bus Voltage Angle ---Load---- -------Generator----- Injected
% No code Mag. Degree MW Mvar MW Mvar Qmin Qmax Mvar
busdata=[1 1 1.04 0.0 0.0 0.0 0.0 0.0 0 0 0
2 2 1.010 0.0 0 -0.8 3.00 88 50 -17 0
3 2 0.985 0.0 40 -1 41 21 60 -10 0
4 0 1.0 0.0 0 0 0 0.0 0 0 0
5 0 1.0 0.0 0 0.0 13 4 0 0 0
6 2 0.98 0.0 0 0.8 75 2 25 -8 0
7 0 1.0 0.0 0 0 0 0.0 0 0 0
8 2 1.005 0.0 450 62.1 150 22 200 -140 0
9 2 0.98 0.0 0 2.2 121 26 9 -3 0
10 0 1.0 0.0 0 0.0 5.0 2 0 0 0
11 0 1.0 0.0 0.0 0.0 0.0 0.0 0 0 0
12 2 1.015 0 310 128.5 377 24 155 -150 0
13 0 1.0 0 0 0 18 2.3 0 0 0
14 0 1 0 0 0 10.5 5.3 0 0 0
15 0 1 0 0 0 22 5 0 0 0
16 0 1 0 0 0 43 3 0 0 0
17 0 1 0 0 0 42 8 0 0 0
18 0 1 0 0 0 27.2 9.8 0 0 0
19 0 1 0 0 0 3.3 0.6 0 0 0
20 0 1 0 0 0 2.3 1 0 0 0
21 0 1 0 0 0 0 0 0 0 0
22 0 1 0 0 0 0 0 0 0 0
23 0 1 0 0 0 6.3 2.1 0 0 0
24 0 1 0 0 0 0 0 0 0 0
25 0 1 0 0 0 6.3 3.2 0 0 0
26 0 1 0 0 0 0 0 0 0 0
27 0 1 0 0 0 9.3 0.5 0 0 0
28 0 1 0 0 0 4.6 2.3 0 0 0
29 0 1 0 0 0 17 2.6 0 0 0
30 0 1 0 0 0 3.6 1.8 0 0 0
31 0 1 0 0 0 5.8 2.9 0 0 0
32 0 1 0 0 0 1.6 0.8 0 0 0
33 0 1 0 0 0 3.8 1.9 0 0 0
34 0 1 0 0 0 0 0 0 0 0
35 0 1 0 0 0 6 3 0 0 0
36 0 1 0 0 0 0 0 0 0 0
37 0 1 0 0 0 0 0 0 0 0
38 0 1 0 0 0 14 7 0 0 0
39 0 1 0 0 0 0 0 0 0 0
40 0 1 0 0 0 0 0 0 0 0
41 0 1 0 0 0 6.3 3 0 0 0
42 0 1 0 0 0 7.1 4.4 0 0 0
43 0 1 0 0 0 2 1 0 0 0
44 0 1 0 0 0 12 1.8 0 0 0
45 0 1 0 0 0 0 0 0 0 0
46 0 1 0 0 0 0 0 0 0 0
47 0 1 0 0 0 29.7 11.6 0 0 0
48 0 1 0 0 0 0 0 0 0 0
49 0 1 0 0 0 18 8.5 0 0 0
50 0 1 0 0 0 21 10.5 0 0 0
51 0 1 0 0 0 18 5.3 0 0 0
52 0 1 0 0 0 4.9 2.2 0 0 0
53 0 1 0 0 0 20 10 0 0 0
54 0 1 0 0 0 4.1 1.4 0 0 0
55 0 1 0 0 0 6.8 3.4 0 0 0
56 0 1 0 0 0 7.6 2.2 0 0 0
57 0 1 0 0 0 6.7 2.0 0 0 0];
% Line code
% Bus bus R X 1/2 B = 1 for lines
% nl nr p.u. p.u. p.u. > 1 or < 1 tr. tap at bus nl
linedata=[1 2 0.0083 0.0280 0.0645 1
2 3 0.0298 0.0850 0.0409 1
3 4 0.0112 0.0366 0.0190 1
4 5 0.0625 0.1320 0.0129 1
4 6 0.0430 0.1480 0.0174 1
6 7 0.0200 0.1020 0.0138 1
6 8 0.0339 0.1730 0.0235 1
8 9 0.0099 0.0505 0.0274 1
9 10 0.0369 0.1679 0.0220 1
9 11 0.0258 0.0848 0.0109 1
9 12 0.0648 0.2950 0.0386 1
9 13 0.0481 0.1580 0.0203 1
13 14 0.0132 0.0434 0.0055 1
13 15 0.0269 0.0869 0.0115 1
1 15 0.0178 0.0910 0.0494 1
1 16 0.0454 0.2060 0.0273 1
1 17 0.0238 0.1080 0.0143 1
3 15 0.0162 0.0530 0.0272 1
4 18 0.0 0.5550 0.0 0.970
4 18 0.0 0.4300 0.0 0.978
5 6 0.0302 0.0641 0.0062 1
7 8 0.0139 0.0712 0.0097 1
10 12 0.0277 0.1262 0.0164 1
11 13 0.0223 0.0732 0.0094 1
12 13 0.0178 0.0580 0.0302 1
12 16 0.0180 0.0813 0.0108 1
12 17 0.0397 0.1790 0.0238 1
14 15 0.0171 0.0547 0.0074 1
18 19 0.4610 0.6850 0.0 1
19 20 0.2830 0.4340 0.0 1
21 20 0.0 0.7767 0.0 1.043
21 22 0.0736 0.1170 0.0 1
22 23 0.0099 0.0152 0.0 1
23 24 0.1660 0.2560 0.0042 1
24 25 0.0 1.1820 0.0 1
24 25 0.0 1.2300 0.0 1
24 26 0.0 0.0473 0.0 1.043
26 27 0.1650 0.2540 0.0 1
27 28 0.0618 0.0954 0.0 1
28 29 0.0418 0.0587 0.0 1
7 29 0.0 0.0648 0.0 0.967
25 30 0.1350 0.2020 0.0 1
30 31 0.3260 0.4970 0.0 1
31 32 0.5070 0.7550 0.0 1
32 33 0.0392 0.0360 0.0 1
34 32 0.0 0.9530 0.0 0.975
34 35 0.0520 0.0780 0.0016 1
35 36 0.0430 0.0537 0.0008 1
36 37 0.0290 0.0366 0.0 1
37 38 0.0651 0.1009 0.0010 1
37 39 0.0239 0.0379 0.0 1
36 40 0.0300 0.0466 0.0 1
22 38 0.0192 0.0295 0.0 1
11 41 0.0 0.7490 0.0 0.955
41 42 0.2070 0.3520 0.0 1
41 43 0.0 0.4120 0.0 1
38 44 0.0289 0.0585 0.0010 1
15 45 0.0 0.1042 0.0 0.955
14 46 0.0 0.0735 0.0 0.900
46 47 0.0230 0.0680 0.0016 1
47 48 0.0182 0.0233 0.0 1
48 49 0.0834 0.1290 0.0024 1
49 50 0.0801 0.1280 0.0 1
50 51 0.1386 0.2200 0.0 1
10 51 0.0 0.0712 0.0 0.930
13 49 0.0 0.1910 0.0 0.895
29 52 0.1442 0.1870 0.0 1
52 53 0.0762 0.0984 0.0 1
53 54 0.1878 0.2320 0.0 1
54 55 0.1732 0.2265 0.0 1
11 43 0.0 0.1530 0.0 0.958
44 45 0.0624 0.1242 0.0020 1
40 56 0.0 1.1950 0.0 0.958
56 41 0.5530 0.5490 0.0 1
56 42 0.2125 0.3540 0.0 1
39 57 0.0 1.3550 0.0 0.980
57 56 0.1740 0.2600 0.0 1
38 49 0.1150 0.1770 0.0015 1
38 48 0.0312 0.0482 0.0 1
9 55 0.0 0.1205 0.0 0.940];
LFYBUS
LFNEWTON
BUSOUT
% This program obtains th Bus Admittance Matrix for power flow solution
% Copyright (c) 1998 by H. Saadat
j=sqrt(-1); i = sqrt(-1);
nl = linedata(:,1); nr = linedata(:,2); R = linedata(:,3);
X = linedata(:,4); Bc = j*linedata(:,5); a = linedata(:, 6);
nbr=length(linedata(:,1)); nbus = max(max(nl), max(nr));
Z = R + j*X; y= ones(nbr,1)./Z; %branch admittance
for n = 1:nbr
if a(n) <= 0 a(n) = 1; else end
Ybus=zeros(nbus,nbus); % initialize Ybus to zero
% formation of the off diagonal elements
for k=1:nbr;
Ybus(nl(k),nr(k))=Ybus(nl(k),nr(k))-y(k)/a(k);
Ybus(nr(k),nl(k))=Ybus(nl(k),nr(k));
end
end
% formation of the diagonal elements
for n=1:nbus
for k=1:nbr
if nl(k)==n
Ybus(n,n) = Ybus(n,n)+y(k)/(a(k)^2) + Bc(k);
elseif nr(k)==n
Ybus(n,n) = Ybus(n,n)+y(k) +Bc(k);
else, end
end
end
clear Pgg
% Power flow solution by Newton-Raphson method
% Copyright (c) 1998 by H. Saadat
ns=0; ng=0; Vm=0; delta=0; yload=0; deltad=0;
nbus = length(busdata(:,1));
for k=1:nbus
n=busdata(k,1);
kb(n)=busdata(k,2); Vm(n)=busdata(k,3); delta(n)=busdata(k, 4);
Pd(n)=busdata(k,5); Qd(n)=busdata(k,6); Pg(n)=busdata(k,7); Qg(n) = busdata(k,8);
Qmin(n)=busdata(k, 9); Qmax(n)=busdata(k, 10);
Qsh(n)=busdata(k, 11);
if Vm(n) <= 0 Vm(n) = 1.0; V(n) = 1 + j*0;
else delta(n) = pi/180*delta(n);
V(n) = Vm(n)*(cos(delta(n)) + j*sin(delta(n)));
P(n)=(Pg(n)-Pd(n))/basemva;
Q(n)=(Qg(n)-Qd(n)+ Qsh(n))/basemva;
S(n) = P(n) + j*Q(n);
end
end
for k=1:nbus
if kb(k) == 1, ns = ns+1; else, end
if kb(k) == 2 ng = ng+1; else, end
ngs(k) = ng;
nss(k) = ns;
end
Ym=abs(Ybus); t = angle(Ybus);
m=2*nbus-ng-2*ns;
maxerror = 1; converge=1;
iter = 0;
% Start of iterations
clear A DC J DX
while maxerror >= accuracy & iter <= maxiter % Test for max. power mismatch
for i=1:m
for k=1:m
A(i,k)=0; %Initializing Jacobian matrix
end, end
iter = iter+1;
for n=1:nbus
nn=n-nss(n);
lm=nbus+n-ngs(n)-nss(n)-ns;
J11=0; J22=0; J33=0; J44=0;
for i=1:nbr
if nl(i) == n | nr(i) == n
if nl(i) == n, l = nr(i); end
if nr(i) == n, l = nl(i); end
J11=J11+ Vm(n)*Vm(l)*Ym(n,l)*sin(t(n,l)- delta(n) + delta(l));
J33=J33+ Vm(n)*Vm(l)*Ym(n,l)*cos(t(n,l)- delta(n) + delta(l));
if kb(n)~=1
J22=J22+ Vm(l)*Ym(n,l)*cos(t(n,l)- delta(n) + delta(l));
J44=J44+ Vm(l)*Ym(n,l)*sin(t(n,l)- delta(n) + delta(l));
else, end
if kb(n) ~= 1 & kb(l) ~=1
lk = nbus+l-ngs(l)-nss(l)-ns;
ll = l -nss(l);
% off diagonalelements of J1
A(nn, ll) =-Vm(n)*Vm(l)*Ym(n,l)*sin(t(n,l)- delta(n) + delta(l));
if kb(l) == 0 % off diagonal elements of J2
A(nn, lk) =Vm(n)*Ym(n,l)*cos(t(n,l)- delta(n) + delta(l));end
if kb(n) == 0 % off diagonal elements of J3
A(lm, ll) =-Vm(n)*Vm(l)*Ym(n,l)*cos(t(n,l)- delta(n)+delta(l)); end
if kb(n) == 0 & kb(l) == 0 % off diagonal elements of J4
A(lm, lk) =-Vm(n)*Ym(n,l)*sin(t(n,l)- delta(n) + delta(l));end
else end
else , end
end
Pk = Vm(n)^2*Ym(n,n)*cos(t(n,n))+J33;
Qk = -Vm(n)^2*Ym(n,n)*sin(t(n,n))-J11;
if kb(n) == 1 P(n)=Pk; Q(n) = Qk; end % Swing bus P
if kb(n) == 2 Q(n)=Qk;
if Qmax(n) ~= 0
Qgc = Q(n)*basemva + Qd(n) - Qsh(n);
if iter <= 7 % Between the 2th & 6th iterations
if iter > 2 % the Mvar of generator buses are
if Qgc < Qmin(n), % tested. If not within limits Vm(n)
Vm(n) = Vm(n) + 0.01; % is changed in steps of 0.01 pu to
elseif Qgc > Qmax(n), % bring the generator Mvar within
Vm(n) = Vm(n) - 0.01;end % the specified limits.
else, end
else,end
else,end
end
if kb(n) ~= 1
A(nn,nn) = J11; %diagonal elements of J1
DC(nn) = P(n)-Pk;
end
if kb(n) == 0
A(nn,lm) = 2*Vm(n)*Ym(n,n)*cos(t(n,n))+J22; %diagonal elements of J2
A(lm,nn)= J33; %diagonal elements of J3
A(lm,lm) =-2*Vm(n)*Ym(n,n)*sin(t(n,n))-J44; %diagonal of elements of J4
DC(lm) = Q(n)-Qk;
end
end
%DX=A\DC';
DX=(inv(A))*DC';
for n=1:nbus
nn=n-nss(n);
lm=nbus+n-ngs(n)-nss(n)-ns;
if kb(n) ~= 1
delta(n) = delta(n)+DX(nn); end
if kb(n) == 0
Vm(n)=Vm(n)+DX(lm); end
end
maxerror=max(abs(DC));
if iter == maxiter & maxerror > accuracy
fprintf('\nWARNING: Iterative solution did not converged after ')
fprintf('%g', iter), fprintf(' iterations.\n\n')
fprintf('Press Enter to terminate the iterations and print the results \n')
converge = 0; pause, else, end
end
if converge ~= 1
tech= (' ITERATIVE SOLUTION DID NOT CONVERGE'); else,
tech=(' Power Flow Solution by Newton-Raphson Method');
end
V = Vm.*cos(delta)+j*Vm.*sin(delta);
deltad=180/pi*delta;
i=sqrt(-1);
k=0;
for n = 1:nbus
if kb(n) == 1
k=k+1;
S(n)= P(n)+j*Q(n);
Pg(n) = P(n)*basemva + Pd(n);
Qg(n) = Q(n)*basemva + Qd(n) - Qsh(n);
Pgg(k)=Pg(n);
Qgg(k)=Qg(n); %june 97
elseif kb(n) ==2
k=k+1;
S(n)=P(n)+j*Q(n);
Qg(n) = Q(n)*basemva + Qd(n) - Qsh(n);
Pgg(k)=Pg(n);
Qgg(k)=Qg(n); % June 1997
end
yload(n) = (Pd(n)- j*Qd(n)+j*Qsh(n))/(basemva*Vm(n)^2);
end
busdata(:,3)=Vm'; busdata(:,4)=deltad';
Pgt = sum(Pg); Qgt = sum(Qg); Pdt = sum(Pd); Qdt = sum(Qd); Qsht = sum(Qsh);
%clear A DC DX J11 J22 J33 J44 Qk delta lk ll lm
%clear A DC DX J11 J22 J33 Qk delta lk ll lm
% This program prints the power flow solution in a tabulated form
% on the screen.
%
% Copyright (C) 1998 by H. Saadat.
%clc
disp(tech)
fprintf(' Maximum Power Mismatch = %g \n', maxerror)
fprintf(' No. of Iterations = %g \n\n', iter)
head =[' Bus Voltage Angle ------Load------ ---Generation--- Injected'
' No. Mag. Degree MW Mvar MW Mvar Mvar '
' '];
disp(head)
for n=1:nbus
fprintf(' %5g', n), fprintf(' %7.3f', Vm(n)),
fprintf(' %8.3f', deltad(n)), fprintf(' %9.3f', Pd(n)),
fprintf(' %9.3f', Qd(n)), fprintf(' %9.3f', Pg(n)),
fprintf(' %9.3f ', Qg(n)), fprintf(' %8.3f\n', Qsh(n))
end
fprintf(' \n'), fprintf(' Total ')
fprintf(' %9.3f', Pdt), fprintf(' %9.3f', Qdt),
fprintf(' %9.3f', Pgt), fprintf(' %9.3f', Qgt), fprintf(' %9.3f\n\n', Qsht)

답변 (4개)

Ahmed Sami
Ahmed Sami 2018년 1월 29일
i'm also interested to know the reason i have the same problem with ieee57 and ieee118.
  댓글 수: 3
AHMED SISTA SOFIANE
AHMED SISTA SOFIANE 2018년 9월 27일
Hi, i have this problem with 39 Bus did you find the answer how to solved Thank you
Iseoluwa Oyedotun
Iseoluwa Oyedotun 2021년 10월 8일
I'm having a similar issue with the IEEE 30 bus not coverging after 100 iterations

댓글을 달려면 로그인하십시오.


Ahmed Sami
Ahmed Sami 2019년 11월 10일
did you find the answer i have the same problem with ieee 14 bus. how to solve? Thank you

kukuh widarsono
kukuh widarsono 2022년 3월 24일
maybe you can increase the generator capacity...

salah Mokred
salah Mokred 2022년 6월 18일
편집: salah Mokred 2023년 3월 10일
This is happened due to the following reasons:
1- You add the loads in the columns of generators and the generators in the columns of load.
2- You should increase the maxiter to 10000 and supposed accel = 1.8.
3- There are some of parallel lines in IEEE57 bus.
You should also see if those data are correct or not by comparing the data with references from the SCI articles.
  댓글 수: 1
SURAJ KUMAR
SURAJ KUMAR 2023년 6월 11일
There are some of parallel lines in IEEE57 bus but print branch power flow only 78 line power finding but two parallel line data not coming 4-18 and 24-25? how solve this probleam..

댓글을 달려면 로그인하십시오.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by