Index exceeds matrix dimensions. error

조회 수: 2 (최근 30일)
giri firmansyah
giri firmansyah 2019년 5월 25일
답변: Walter Roberson 2019년 5월 26일
Hello, im trying to find location of distributed generation to reduce line losses.
I kept getting this error.
Index exceeds matrix dimensions.
Error in DG_PSO (line 149)
lk = nbus+l-ngs(l)-nss(l)-ns;
this is the m-file
%PSO Algoritm for DG position assignment
%Loss Reductin
clear all
clc
%====================Random DG position &active power =====================
for i=1:10
x(i,1)=ceil(38*rand);
x(i,2)=6.5*rand;
end
for q=1:100
pop_size=10;
for h=1:pop_size
basemva=10000; accuracy=0.0001; maxiter= 1000;
% Bus Bus Volt. Angle ---Load--- ---------Generator-------- Injected
% no code mag. Deg. KW KVAR KW KVAR Qmin Qmax KVAR
% Bus Bus Volt. Angle ---Load--- ---------Generator-------- Injected
% no code mag. Deg. KW KVAR KW KVAR Qmin Qmax KVAR
busdata=[1 1 1.00000 0.0 0 0 3480.00 2225.3 0 0 0
2 2 0.99997 0.0 0 0 0 0 0 0 0
3 0 0.99965 0.0 114.917 73.229 0 0 0 0 0
4 0 0.98435 0.0 7.043 4.392 0 0 0 0 0
7 0 0.99565 0.0 422.310 270.023 0 0 0 0 0
8 0 0.98513 0.0 54.121 34.381 0 0 0 0 0
12 0 0.99259 0.0 351.813 223.411 0 0 0 0 0
13 0 0.99077 0.0 53.176 33.745 0 0 0 0 0
15 0 0.98910 0.0 107.085 68.041 0 0 0 0 0
21 0 0.98769 0.0 32.355 20.642 0 0 0 0 0
23 0 0.98683 0.0 27.151 17.242 0 0 0 0 0
25 0 0.98594 0.0 24.188 15.320 0 0 0 0 0
27 0 0.98447 0.0 26.398 16.754 0 0 0 0 0
30 0 0.98437 0.0 0 0 0 0 0 0 0
31 0 0.98436 0.0 4.339 2.700 0 0 0 0 0
34 0 0.98430 0.0 14.703 9.233 0 0 0 0 0
36 0 0.98428 0.0 8.430 5.304 0 0 0 0 0
38 0 0.98417 0.0 70.150 44.689 0 0 0 0 0
42 0 0.98407 0.0 34.487 22.049 0 0 0 0 0
44 0 0.98428 0.0 14.345 9.006 0 0 0 0 0
47 0 0.98413 0.0 24.369 15.438 0 0 0 0 0
49 0 0.98405 0.0 28.335 18.015 0 0 0 0 0
56 0 0.99021 0.0 321.030 202.724 0 0 0 0 0
62 0 0.98904 0.0 20.236 12.770 0 0 0 0 0
63 0 0.98779 0.0 91.152 57.598 0 0 0 0 0
68 0 0.98966 0.0 139.110 88.631 0 0 0 0 0
71 0 0.98926 0.0 196.816 125.697 0 0 0 0 0
74 0 0.98886 0.0 68.711 43.912 0 0 0 0 0
76 0 0.98878 0.0 29.306 18.644 0 0 0 0 0
77 0 0.98874 0.0 41.877 26.944 0 0 0 0 0
80 0 0.98913 0.0 44.227 28.152 0 0 0 0 0
83 0 0.98953 0.0 132.245 83.935 0 0 0 0 0
87 0 0.98918 0.0 123.454 77.955 0 0 0 0 0
90 0 0.98720 0.0 149.767 94.922 0 0 0 0 0
92 0 0.98707 0.0 46.873 29.607 0 0 0 0 0
96 0 0.98696 0.0 86.364 54.695 0 0 0 0 0
97 0 0.98629 0.0 475.504 301.828 0 0 0 0 0
101 0 0.98611 0.0 64.017 40.832 0 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.0007252 0.0002611 0 1
1 7 0.0158750 0.0057161 0 1
2 3 0.0028587 0.0003062 0 1
7 12 0.0026255 0.0009454 0 1
8 27 0.0185094 0.0019828 0 1
12 13 0.0101557 0.0010879 0 1
12 56 0.0058843 0.0006303 0 1
12 83 0.0051502 0.0005517 0 1
13 15 0.0017798 0.0001906 0 1
15 62 0.0023470 0.0013226 0 1
15 63 0.0194216 0.0020805 0 1
23 25 0.0073637 0.0007888 0 1
25 8 0.0243159 0.0026048 0 1
27 30 0.0240267 0.0025738 0 1
27 38 0.0252028 0.0027001 0 1
27 44 0.0124027 0.0013286 0 1
30 4 0.0043781 0.0018003 0 1
30 31 0.0014747 0.0006064 0 1
30 34 0.0124983 0.0051392 0 1
34 36 0.0020001 0.0008224 0 1
38 42 0.0086271 0.0035474 0 1
44 47 0.0213459 0.0022866 0 1
47 49 0.0176874 0.0072730 0 1
56 68 0.0052061 0.0018746 0 1
63 21 0.0164405 0.0017611 0 1
63 23 0.0167519 0.0017945 0 1
68 71 0.0060980 0.0021957 0 1
71 74 0.0030590 0.0003277 0 1
71 80 0.0176084 0.0018862 0 1
74 76 0.0028757 0.0011825 0 1
74 77 0.0024333 0.0010006 0 1
83 87 0.0063515 0.0006804 0 1
83 90 0.0141268 0.0015133 0 1
90 92 0.0198999 0.0021317 0 1
90 96 0.0055284 0.0005922 0 1
90 97 0.0163568 0.0030762 0 1
97 101 0.0062291 0.0006673 0 1];
LFYBUS;
% 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';
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
%===============================power loss================================
SLT = 0;
for n = 1:nbus
busprt = 0;
for L = 1:nbr;
if busprt == 0
busprt = 1;
else, end
if nl(L)==n k = nr(L);
In = (V(n) - a(L)*V(k))*y(L)/a(L)^2 + Bc(L)/a(L)^2*V(n);
Ik = (V(k) - V(n)/a(L))*y(L) + Bc(L)*V(k);
Snk = V(n)*conj(In)*basemva;
Skn = V(k)*conj(Ik)*basemva;
SL = Snk + Skn;
SLT = SLT + SL;
elseif nr(L)==n k = nl(L);
In = (V(n) - V(k)/a(L))*y(L) + Bc(L)*V(n);
Ik = (V(k) - a(L)*V(n))*y(L)/a(L)^2 + Bc(L)/a(L)^2*V(k);
Snk = V(n)*conj(In)*basemva;
Skn = V(k)*conj(Ik)*basemva;
SL = Snk + Skn;
SLT = SLT + SL;
else, end
end
end
SLT = SLT/2;
f(h,1)=real(SLT);
busdata(x(h,1),7)=0;
end
%============================PSO algoritm(main program)====================
if q==1
g_position=f(1); p_position=f; iteration=q; p_best(h,:)=x(h,:); g_best=x(1,:);
velocity=zeros(10,2);
end
w=0.9;
c1=0.7;
c2=0.7;
r1=rand(10,1);
r2=rand(10,1);
for y=1:10
if p_position(y)>f(y)
p_position(y)=f(y);
p_best(y,:)=x(y,:);
end
end
for u=1:10
if g_position>p_position(u)
g_position=p_position(u);
g_best=x(u,:);
iteration=q;
end
end
for k=1:10
velocity(k,:)=w*velocity(k,:)+r1(k,1)*c1.*(p_best(k,:)-x(k,:))+r2(k,1)*c2.*(g_best-x(k,:));
x(k,1)=ceil(x(k,1)+velocity(k,1));
x(k,2)=x(k,2)+velocity(k,2);
if x(k,1)>33|x(k,1)<2
x(k,1)=ceil(2+31*rand);
end
if x(k,2)<0
x(k,2)=2;
end
end
end
%=========================best position of DG==============================
power_loss=g_position;
g_best , g_position
can anybody help me to solve this error?
any help would be appreciated, thank you

답변 (1개)

Walter Roberson
Walter Roberson 2019년 5월 26일
Your code is using the bus numbers, up to 101, as indices into arrays that as sized according to the number of busses. In places you should be dropping down from the bus number to the bus index. You could create a de-reference vector like
bus_numbers = busdata(:,1);
bidx(bus_numbers) = 1:length(bus_numbers);

카테고리

Help CenterFile Exchange에서 Arithmetic Operations에 대해 자세히 알아보기

제품


릴리스

R2017a

Community Treasure Hunt

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

Start Hunting!

Translated by