clear all
Ex30_2013_Rev

MontaY_Rev2

%Tensões pré-falta
%Vpre=ones(nb,1)*j;
Vpre=ones(nb,1);
aa=1*exp(j*120*pi/180);
A=[ 1 1 1; 1 aa*aa aa; 1 aa aa*aa];
Yshunt=bl/2;
y1=1./(j*x1);
y2=1./(j*x2);
y0=1./(j*x0);
delete saida.out
diary saida.out
% Cálculo de Curto Trifásico
% Corrente de Falta
Ifff=sparse(nb,1);
Ifff=Vpre./diag(Z1);

modu=abs(Ifff);
angu=(180/pi)*angle(Ifff);

disp('      Cálculo de Curto Trifásico                      ')
disp('    barra    Corrente de Falta  pu   ')
disp([biv  modu angu ]);

% Tensão pós-falta
for i=1:nb  % falta na barra i
   respline=[];
   resp=[];
   for l=1:nb
   Vpos{i}(l,1)= Vpre(l,1)-Vpre(i,1)*(Z1(l,i)/Z1(i,i));
   x=A*[ 0;Vpos{i}(l,1); 0];
   mod{l}=abs(x);
   ang{l}=(180/pi)*angle(x);
   resp=[ resp; [mod{l}(1,1) ang{l}(1,1) mod{l}(2,1) ang{l}(2,1) mod{l}(3,1) ang{l}(3,1)]];
   end
   
bar_i=Af'*Vpos{i};
bar_j=Afo'*Vpos{i};

Itotalij= y1.*bar_i - y1.*bar_j + bar_i.*Yshunt;
Itotalji= y1.*bar_j -y1.*bar_i + bar_j.*Yshunt;

for line=1:nl
x=A*[ 0;Itotalij(line,1); 0];
modline{line}=abs(x);
angline{line}=(180/pi)*angle(x);
respline=[ respline; [modline{line}(1,1) angline{line}(1,1) modline{line}(2,1) angline{line}(2,1)...
           modline{line}(3,1) angline{line}(3,1)]];
end
   
   disp(['    barra    Tensão Pós-falta para  Falta trifásica na barra  ' num2str(i)   ]  )
   disp(['                  fase a               fase b              fase c '  ]  )
   disp([biv resp])
   
   disp(['    linha   Correntes  Pós-falta para  Falta trifásica na barra  ' num2str(i)   ]  )
   disp(['                  fase a               fase b              fase c '  ]  )
   disp([liv respline])
 
end

diary off
