clc
//////////Variables de entrada//
LARGO=450
G=40*10^-9
f=50
S=100
V=66
seccion=725
RY=0.0417
n=3
r=(sqrt(seccion/%pi))*(0.001)                                               /// Radio de mm2 a m
while (n>4 || n<0)

    n=input('Número no permitido,Ingrese el número de subconductores (1 a 4) ')
end
d=0.4
D12=4.8052
D13=4.8052
D23=9.4

//////Condiciones

if n==1 then reql=0.7788007831*r                                            //// Condicion para REQL
   elseif n==2 then reql=sqrt(0.778800783*r*d)
   elseif n==3 then reql=(0.778800783*r*d^2)^(1/3)
   elseif n==4 then reql=1.09*((0.778800783*r*d^3)^(1/4))
   else reql=r
end

if n==1 then reqc=r*d               //// Condicion para REQC
elseif n==2 then reqc=sqrt(r*d)
elseif n==3 then reqc= (r*d^(2))^(1/3)
elseif n==4 then reqc=1.09*((r*(d^(3)))^(1/4))
else reqc=r
end

//Correccion de temperatura

T=228
TF=35
RO=(((T+TF)*RY)/(T+25))/n

////ecuaciones
W=(2*%pi*f)             ///W
E0=8.8542*10^(-12)          //:)
//ecuaciones de inductancia
Deq=(D12*D13*D23)^(1/3)
IND=2*10^(-7)*log(Deq/reql)*1000
ZL=W*IND

//ecuaciones de capacitancia
CAP=((2*%pi*E0)/log(Deq/reqc))*1000
YC=CAP*W

//ecuaciones para matriz exacta
Z=complex(RO,ZL)
Y=complex(G,YC)
gama=((complex(RO,W*IND))*(complex(G,W*CAP)))^(0.5)
ZC=((complex(RO,W*IND))/(complex(G,W*CAP)))^(0.5)

//
ZB=(V^2/S)

//Ecuaciones para la matriz PI //
AE=cosh(gama*LARGO)
BE=ZC*(sinh(gama*LARGO))
CE=(sinh(gama*LARGO))/ZC
DE=AE
BPU=BE/ZB
CPU=CE*ZB
MT=[AE,BE;CE,DE]
MTPU=[AE,BPU;CPU,DE]

//Ecuaciones para la matriz CORTA PU//
ZA=(complex(RO,ZL)*LARGO)
YA=(complex(0,YC)*LARGO)
AC=((ZA*YA)/2)+1
BC=ZA/ZB
CC=YA*(1+((ZA*YA)/4))*ZB
DC=AC
MC=[AC,BC;CC,DC]
BC2=BC/2
CC2=CC*2
MC2=[AC,BC2;CC2,DC]
/*
disp(ZA, 'Z aproximada')
disp(YA, 'Y aproximada')
disp(RO, 'resistencia corregida')


////IMPRESIONES

disp(Deq, 'Distancia equivalente')
disp(reqc, 'reqc')
disp(reql, 'reql')
disp(CAP, 'Capacitancia')
disp(IND, 'Inductancia')
disp(RO, 'Resistencia arreglada')
disp(gama, 'gama')
disp(ZC, 'zc')
disp(MT, 'Matriz')
disp(MTPU, 'MTPU')
disp(MC, 'Matriz corta') */
CR=2
BPU2=(BE/ZB)/CR
CPU2=CE*ZB*CR
MT2PU=[AE,BPU2;CPU2,DE]
V2=complex(1,0)
I2=complex(0.96,-0.27999)



//rendimiento sin falla
V1=(AE*V2)+(BPU2*I2)
I1=(CPU2*V2)+(DE*I2)

aI2=real(I2)
bI2=imag(I2)
aI1=real(I1)
bI1=imag(I1)
rendimiento=real(V2*complex(aI2,-bI2))/real(V1*complex(aI1,-bI1))*100

//Rendimiento con falla
V12=(AE*V2)+(BPU*I2)
I12=(CPU*V2)+(DE*I2)
aI12=real(I12)
bI12=imag(I12)
rendimiento2=real(V2*complex(aI2,-bI2))/real(V12*complex(aI12,-bI12))*100

//ferranti
FERR=1/(sqrt(real(AE)^2+imag(AE)^2))  //SIN  FALLA
FERR2=1/(sqrt(real(AC)^2+imag(AC)^2)) //CON FALLA

//regulacion
reg=((sqrt(real(V1)^2+imag(V1)^2))-(sqrt(real(V2)^2+imag(V2)^2)))/(sqrt(real(V1)^2+imag(V1)^2)) //SIN FALLA
reg2=((sqrt(real(V12)^2+imag(V12)^2))-(sqrt(real(V2)^2+imag(V2)^2)))/(sqrt(real(V12)^2+imag(V12)^2)) // CON FALLA

//Potencia aparente
S1=V1*(complex(real(I1),-imag(I1))) // sin falla
S2=V12*(complex(real(I12),-imag(I12))) // sin falla



/// impresiones
disp(MT2PU, 'Matriz Exacta')
disp(MC2, 'Matriz Aproximada')
disp(FERR, 'Efecto Ferranti matriz exacta')
disp(FERR2, 'Efecto Ferranti matriz aproximada')
disp(rendimiento, 'Rendimiento Matriz exacta')
disp(reg*100, 'Regulacion matriz exacta')
disp(rendimiento2, 'Rendimiento matriz exacta con falla')
disp(reg2*100, 'Regulacion matriz exacta con Falla')
disp(S1, 'Potencia Aparente sin falla')
disp(S2, 'Potencia Aparente Con falla')
disp(V1, 'Voltaje 1')
disp(V2, 'Voltaje 2')
disp(I1, 'Corriente 1')
disp(I2, 'Corriente 2')


/*
disp(MTPU, 'Matriz con falla')


disp(rendimiento2, 'Rendimiento con falla')
disp(V1, 'Voltaje emisor')
disp(I1, 'CORRIENTE EMISOR')
disp(V12, 'Voltaje emisor con falla')
disp(I12, 'CORRIENTE EMISOR con falla')

*/