error:subscript indicates must be real or positive...

HI my problem is when I run this program, in first loop(*for ui=1:2) there is'nt any error , but in second loop this error appears in line :

****vb_phasor=[abs(vb) angle(vb)*180/pi]
   'subscript must be ....'

I was wondering if this program is problem , why it runs in first loop?

pp.this code related to load flow for radial network when we have DG .

clear all;
format short
G;
tic
m=load('loaddata14bus.m');
l=load('linedata14bus.m');
o=load('DGcapacity.m');
br=length(l);
no=length(m);
MVAb=100;
KVb=20;
Zb=(KVb^2)/MVAb;

% Per unit Values
for i=1:br
   R(i,1)=(l(i,4));
   X(i,1)=(l(i,5));
end
ui=1;

while  ui<=2     ****%%THIS LINE

for i=1:no
    P(i,1)=((m(i,2)-o(i,ui))/(1000*MVAb));
    Q(i,1)=((m(i,3)-m(i,5))/(1000*MVAb));
end
R;
X;
P
Q
C=zeros(br,no);

for i=1:br
    a=l(i,2);
    b=l(i,3);
    for j=1:no
        if a==j
            C(i,j)=-1;
        end
        if b==j
            C(i,j)=1;
        end
    end
end
C;
e=1;
for i=1:no
    d=0;
    for j=1:br
        if C(j,i)==-1
            d=1;
        end
    end
    if d==0
        endnode(e,1)=i;
        e=e+1;
    end
end
endnode
h=length(endnode);
for j=1:h
    e=2;

      f=endnode(j,1);
     % while (f~=1)
     for s=1:no
       if (f~=1)
         k=1;  
         for i=1:br
             if ((C(i,f)==1)&&(k==1))
                  f=i;
                  k=2;
             end
         end
         k=1;
         for i=1:no
             if ((C(f,i)==-1)&&(k==1));
                  f=i;
                  g(j,e)=i;
                  e=e+1;
                  k=3;
             end            
         end
       end
     end
  end
  for i=1:h
      g(i,1)=endnode(i,1);
  end
  g;
  w=length(g(1,:));
  for i=1:h
      j=1;
      for k=1:no 
          for t=1:w
              if g(i,t)==k
                  g(i,t)=g(i,j);
                  g(i,j)=k;
                  j=j+1;
               end
           end
      end
  end
  g;
  for k=1:br
      e=1;
      for i=1:h
          for j=1:w-1
              if (g(i,j)==k) 
                  if g(i,j+1)~=0
                      adjb(k,e)=g(i,j+1);            
                      e=e+1;
                  else
                      adjb(k,1)=0;
                  end
               end
          end
      end
  end
  adjb;
  for i=1:br-1
      for j=h:-1:1
          for k=j:-1:2
              if adjb(i,j)==adjb(i,k-1)
                  adjb(i,j)=0;
              end
          end
      end
  end
  adjb;
  x=length(adjb(:,1));
  ab=length(adjb(1,:));
  for i=1:x
      for j=1:ab
          if adjb(i,j)==0 && j~=ab
              if adjb(i,j+1)~=0
                  adjb(i,j)=adjb(i,j+1);
                  adjb(i,j+1)=0;
              end
          end
          if adjb(i,j)~=0
              adjb(i,j)=adjb(i,j)-1;
          end
      end
  end
  adjb;
  for i=1:x-1
      for j=1:ab
          adjcb(i,j)=adjb(i+1,j);
      end
  end
  b=length(adjcb);

% voltage current program

for i=1:no
    vb(i,1)=1;
end

for s=1:10
for i=1:no
    nlc(i,1)=conj(complex(P(i,1),Q(i,1)))/(vb(i,1));
end

nlc;

for i=1:br
    Ibr(i,1)=nlc(i+1,1);
end

Ibr;

xy=length(adjcb(1,:));

for i=br-1:-1:1

      for k=1:xy

          if adjcb(i,k)~=0
              u=adjcb(i,k);
              %Ibr(i,1)=nlc(i+1,1)+Ibr(k,1);
              Ibr(i,1)=Ibr(i,1)+Ibr(u,1);
          end
      end      
  end

Ibr;

for i=2:no
      g=0;
      for a=1:b 
          if xy>1
            if adjcb(a,2)==i-1 
                u=adjcb(a,1);
                vb(i,1)=((vb(u,1))-((Ibr(i-1,1))*(complex((R(i-1,1)),X(i-1,1)))));
                g=1;
            end
            if adjcb(a,3)==i-1 
                u=adjcb(a,1);
                vb(i,1)=((vb(u,1))-((Ibr(i-1,1))*(complex((R(i-1,1)),X(i-1,1)))));
                g=1;
            end
          end
        end
        if g==0
            vb(i,1)=((vb(i-1,1))-((Ibr(i-1,1))*(complex((R(i-1,1)),X(i-1,1)))));
        end
end

s=s+1;
end
nlc;

Ibr;

vb

vb_phasor=[abs(vb) angle(vb)*180/pi]                *****%%THIS LINE HAVE ERROR

toc;

for i=1:no

      va(i,2:3)=vb_phasor(i,1:2);
  end
  for i=1:no
      va(i,1)=i;
  end
  va;

Ibr_phase=[abs(Ibr) angle(Ibr)*180/pi];
PLoss_tot(1,1)=0;
QLoss_tot(1,1)=0;

% losses
for f=1:br
    Pl(f,1)=(Ibr_phase(f,1)^2)*R(f,1);
    Ql(f,1)=X(f,1)*(Ibr_phase(f,1)^2);
    PLoss_tot(1,1)=PLoss_tot(1,1)+Pl(f,1);
    QLoss_tot(1,1)=QLoss_tot(1,1)+Ql(f,1);
end

Plosskw=(Pl)*100000
Qlosskw=(Ql)*100000
PLoss_tot=(PLoss_tot)*100000
QLoss_tot=(QLoss_tot)*100000

voltage = vb_phasor(:,1)
angle = vb_phasor(:,2)*(pi/180)
xlswrite('Result',angle,'PowerFlow','G3:G16')
xlswrite('Result',voltage,'PowerFlow','F3:F16')
xlswrite('Result',Plosskw,'PowerFlow','D3:D16')
xlswrite('Result',Qlosskw,'PowerFlow','E3:E16')
xlswrite('Result',P,'PowerFlow','B3:B16')
xlswrite('Result',Q,'PowerFlow','C3:C16')
ui=ui+1;

end

[per isakson]

Without input data we cannot debug this code.

What is the purpose of lines like

    R;
    X;

You need to use the debugging tools of Matlab. Here are some links:

• Debugging in MATLAB
• Debugging session
• Debugging Approaches

[adele]

that lines related to loadflow.resistance & reactance of line, P:injection power to each buses...

Answer 1

[Star Strider]

This line is the problem:

angle = vb_phasor(:,2)*(pi/180)

You use ‘angle’ as a variable (this is called ‘overshadowing’), so the second iteration of the loop, MATLAB considers ‘angle’ to be a variable and not the built-in function, angle. Rename ‘angle’ to something like ‘vb_angle’ (or something else that makes sense in the context of your code). I only see it otherwise in the xlswrite statement just below it, so an extensive revision won’t be necessary. (You might also consider renaming ‘voltage’ as ‘vb_voltage’ for consistency, but it’s not causing any problem with your code.)

[adele]

I did what u commented me,but this error appears: Undefined function 'vb_angle' for input arguments of type 'double'. Error in check1 (line 203) vb_phasor=[abs(vb) vb_angle(vb)*180/pi]

[Star Strider]

I apologise for not making my Answer clearer.

That line should be:

vb_phasor = [abs(vb) angle(vb)*180/pi]

[adele]

I want to thank u...it's working

[Star Strider]

My pleasure!

Answer 2

[Image Analyst]

See the FAQ: http://matlab.wikia.com/wiki/FAQ#How_do_I_fix_the_error_.22Subscript_indices_must_either_be_real_positive_integers_or_logicals..22.3F

If that doesn't explain it, then step through the code and figure out what values your indexes are. You can do that just as well as I can. If not, see this