Index exceeds matrix dimensions. error

9 次查看（过去 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 Center 和 File Exchange 中查找有关 Parallel for-Loops (parfor) 的更多信息

标签

产品


版本

R2017a

Community Treasure Hunt

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

Start Hunting!

Translated by