whay im getting this erro

clc;
clear all;

tic
Egenerators=3;   
% 1.unit no  2.Pmin  3.Pmax  4.SU  5.SD  6.MUT  7.MDT  8.U0  9.S0  10.RU  11.RD  12.UNIT   13.a($/h)  14.b($/MWh)  15.c($/(MW)2h)

ps.G=[
1  20  175  30  30  3  3  0  1  50  50  1  100  2.45  0.0012;
2  40  300  50  60  3  3  0  1  60  60  2  120  2.32  0.0010;
3  50  500  60  50  4  4  3  0  80  80  3  150 2.10  0.0015;
];

EHgenerator=2;
% 1.unit NO  2.Vmin  3.Vmax  4.VO  5.VT  6.Qmin  7.Qmax  8.Pmin  9.Pmax  10.C1  11.C2  12.C3  13.C4  14.C5  15.C6
ps.H=[
 1  100  240  170  170  10  30  0  500  -0.0016  -0.3  0.014  0.55  5.5  -40;
 2  70  160  120  140  6  20  0  500  -0.003  -0.31  0.027  1.44  14  -90;
];

Hiflow=2;
% % reservoir inflow from hour 1 to 24 
% [8.1  8.2  4  2  3  4  3  2  1  1  1  2  4  3  3  2  2  2  1  1  2  2  1  0   ;
% 2.8  2.4  1.6  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0   ;
% ];

%7条线路
Ebranches=9;
% 1. line.no    2.from bus    3.to bus     4. X (p.u)    5. transmission limit    

ps.Ebranch=[
1  1  2  0.12  40;
2  1  8  0.38  30;
3  2  3  0.24  45;
4  2  6  0.08  20;
5  3  4  0.4     30;
6  4  5  0.08  60;
7  5  6  0.24  70;
8  6  7  0.16  50;
9  7  8  0.08  50;
];

%6个节点
Ebuses=8 ;
%1. Node number 2-25: load for each node at 1-24
ps.Ebus =[     
1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  ;
2  101.256  96.156  75.576  59.1  60.072  61.56  85.404  112.548  118.2  129.768  131.616  121.74  119.196  113.676  131.628  138.06  124.728  130.38  138.276  143.544  144.048  130.308  125.34  117.048  ;
3  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0   ;
4  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  ;
5  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  ;
6  243  230.784  181.38  141.852  144.168  147.744  204.96  270.12  283.692  311.436  315.888  292.176  286.08  272.808  315.9  331.344  299.352  312.912  331.86  344.496  345.72  312.732  300.816  280.908  ;
7   60.756  57.696  45.348  35.46  36.036  36.936  51.24  67.536  70.92  77.856  78.972  73.044  71.52  68.208  78.972  82.836  74.832  78.228  82.968  86.124  86.424  78.18  75.204  70.224  ;
8   0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0   ;
];

Ewinds=2;
% wind power 24 hours output MW
ps.PWforecast=[193.9  195.3  197.1  189.66  200.27  189.34  163.28  156.7  148.7  149.5  143.44  128.6  115.96  110.48  113.22  118.16  122.02  128.6  141.8  152.24  158.18  168.62  178.52  185.62;
    193.9  195.3  197.1  189.66  200.27  189.34  163.28  156.7  148.7  149.5  143.44  128.6  115.96  110.48  113.22  118.16  122.02  128.6  141.8  152.24  158.18  168.62  178.52  185.62;
  ];

%%%forms a susceptance matrix
B0=zeros(Ebuses);

% first calculate self-admittance 
for i=1:Ebuses
  for k=1:Ebranches    
      if ps.Ebranch(k,2)==i||ps.Ebranch(k,3)==i
          B0(i,i)=B0(i,i)+1/ps.Ebranch(k,4); 
      end
  end
end

% %再算互导纳% reciprocal
for k=1:Ebranches
    for i=1:Ebuses
        for j=1:Ebuses
           if (ps.Ebranch(k,2)==i&&ps.Ebranch(k,3)==j)&&i~=j
                B0(i,j)=B0(i,j)-1/ps.Ebranch(k,4);
           end
        end
    end
end

% form a symmetric admittance matrix
for i=1:Ebuses    %矩阵左下半部分
    for j=1:i
        if i~=j&&B0(i,j)~=0
            B0(j,i)=B0(i,j);
        end
    end
end
for i=1:Ebuses   % upper right part of the matrix
    for j=i:Ebuses
        if i~=j&&B0(i,j)~=0
            B0(j,i)=B0(i,j);
        end
    end
end

%%电抗矩阵%%reactance matrix
X=zeros(Ebuses);
X=inv(sym(B0));
% X=inv(B0);

% % power generation transfer factor  I DON'T HAVE IT i skip it
GSDF=zeros(Ebranches,Ebuses); 

for i=1:Ebranches
    for j=1:Ebuses
%         GSDF(i,j)=(X(ps.Ebranch(i,2),j)-X(ps.Ebranch(i,3),j))/ps.Ebranch(i,4);
          GSDF(i,j)=(X(ps.Ebranch(i,2),j)-X(ps.Ebranch(i,3),j))/ps.Ebranch(i,3);
    end
end 

Error using symengine (line 59)
Cannot define a matrix over 'Dom::ExpressionField()'.

Error in sym/subsref (line 773)
        B = mupadmex('symobj::subsref',L.s,privformat(R_tilde));

Error in Untitled (line 114)
          GSDF(i,j)=(X(ps.Ebranch(i,2),j)-X(ps.Ebranch(i,3),j))/ps.Ebranch(i,3);

Answer 1

[Walter Roberson]

Your B0 becomes an 8 x 8 matrix with rank 7. You then inv(sym(B0)) which returns FAIL because of the singularity. Then later when you try to index the FAIL that results from that, that is giving you the error you see.

The work around is to ensure that your matrix is not singular.

[Eltalib Mohamed]

how can i solve it please

[Walter Roberson]

You need to recheck your information about which elements are connected to which other elements.

It might help to look at https://www.ikbooks.com/home/samplechapter?filename=68_Sample_Chapter.pdf