Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad...

25
Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión de Energía IV. Enero-Marzo 2013. Prof. José Manuel Aller. Alumno: José M. Macías D. 09-10466. Tarea N° 2. Una máquina sincrónica de rotor liso de 100 MVA de potencia nominal, 10 kV, fp nominal 0.85, un par de polos, 60 Hz, corriente de campo nominal 300 A, tiene una reactancia de cortocircuito de 1,0 pu. La reactancia de dispersión es de 0.2 pu. 1.- Calcule la máxima potencia reactiva que puede entregar la máquina como condensador sincrónico. A través del siguiente código se obtiene la característica de saturación de la máquina: function plsaturation(Lm0, Lmsat, PsiT, fT) Lm0=2; Lmsat=0.2; PsiT=0.93; fT=1; iT = 1/Lm0*PsiT; Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; Mi = (1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); im = (Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) - PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) - log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; plot(im,Psim); ylabel('\Psi_m'); xlabel('i_m'); grid on; Obteniendo: 0 0.5 1 1.5 2 2.5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 m i m

Transcript of Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad...

Page 1: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

CT-4311 Conversión de Energía IV.

Enero-Marzo 2013.

Prof. José Manuel Aller.

Alumno: José M. Macías D. 09-10466.

Tarea N° 2.

Una máquina sincrónica de rotor liso de 100 MVA de potencia nominal, 10 kV, fp

nominal 0.85, un par de polos, 60 Hz, corriente de campo nominal 300 A, tiene una

reactancia de cortocircuito de 1,0 pu. La reactancia de dispersión es de 0.2 pu.

1.- Calcule la máxima potencia reactiva que puede entregar la máquina como

condensador sincrónico.

A través del siguiente código se obtiene la característica de saturación de la

máquina: function plsaturation(Lm0, Lmsat, PsiT, fT) Lm0=2; Lmsat=0.2; PsiT=0.93; fT=1; iT = 1/Lm0*PsiT; Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; Mi = (1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); im = (Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; plot(im,Psim); ylabel('\Psi_m'); xlabel('i_m'); grid on;

Obteniendo:

0 0.5 1 1.5 2 2.50

0.2

0.4

0.6

0.8

1

1.2

1.4

m

im

Page 2: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Se hallará una recta que se ajuste a la zona lineal de la característica de

magnetización de la máquina:

En el punto nominal tenemos:

Teniendo estos resultados podemos obtener todos los demás datos del punto nominal:

Con Dn queda definido el eje en cuadratura y por consiguiente el eje directo, ahora

hallaremos las proyecciones del vector Een sobre dichos ejes para saber el grado de

saturación de cada uno de ellos:

0 0.5 1 1.5 2 2.50

0.2

0.4

0.6

0.8

1

1.2

1.4

m

im

Característica Saturada

Característica Lineal

Page 3: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Nos resulta en la zona lineal, por lo tanto no se satura la Xq, ahora verificaremos el

Eq para saber si debemos modificar Xd:

Resulta por encima de la zona lineal, así que hallaremos el factor de saturación para

poder corregir el valor de Xd:

Sabemos que en el punto nominal la Ef es la máxima, por lo que la corriente de

campo será máxima de igual forma:

Donde Id será la proyección de la corriente del estator sobre el eje directo:

Page 4: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Quedando la fuerza electromotriz máxima como:

Ahora para calcular la potencia reactiva máxima para la condición de condensador

sincrónico, empezaré asumiendo que la corriente de campo será la máxima:

El Ed es igual a cero, por lo tanto no hace falta calcular la saturación del eje en

cuadratura. Ahora calcularemos el grado de saturación del eje directo para verificar si es el

mismo que utilizamos para este punto de operación (Condensador sincrónico):

Quedando el nuevo Ef como:

Utilizando una rutina en MATLAB se obtendrá el valor en el que se estabiliza Sd

para poder hallar el valor de la potencia reactiva máxima que suministra la máquina en la

condición de operación de condensador sincrónico:

Page 5: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Código de MATLAB 1 hold off clear all clc Lm0=2; Lmsat=0.2; PsiT=.93; fT=1; iT = 1/Lm0*PsiT; Psimax =

4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; Efmax=1.3889; Ven=1; delta=0;s=2.0532; ss=0; Xd=1; Xdisp=0.2; Xq=1; n=-1; hold on while(abs(s-ss)>.0001) n = n+1; s = (s+ss)/2; Xds = Xd/s+((s-1)/s)*Xdisp; Q = (Efmax/Xds)-(1/Xds); Ie = -j*Q; %D = Ven + j*Xq*Ie; Ee = Ven + j*Xdisp*Ie; Eq = Ee; ifl = abs(Eq)/1.8; Psim = abs(Eq); Mi = (1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifs = (Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+.5/tauT*(log(1+(tauT*PsiT)^2) - log(1+tauT^2*(Psim-

PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ss = ifs/ifl; if ss<1 ss=1;

5 10 15 20 25 30 35

0.5

1

1.5

2

2.5

Número de iteraciones

Sd (

azul) E

f (v

erd

e)

Page 6: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía end Efmax = 1.8*1.5843/ss; scatter(n,s,'b','d','filled') scatter(n,Efmax,'g','d','filled') grid on end

2.- La corriente de campo máxima.

Como ya se determinó en la pregunta anterior a partir del punto nominal, la

corriente de campo máxima será:

3.- La corriente de campo mínima para potencia activa nominal.

Para este caso supondremos tensión de estator nominal (Ve=1 p.u.) y variaremos el

valor de la fuerza electromotriz del campo y el ángulo de carga (δ) hasta conseguir el

mínimo valor de corriente de campo para el que podamos suministrar la potencia activa

nominal

Para se logra suministrar la potencia activa

nominal:

Page 7: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Quedando la corriente de campo mínima para suministrar la potencia activa nominal

en:

Código de MATLAB 2 hold off clear all clc Lm0=2; Lmsat=0.2; PsiT=.93; fT=1; iT = 1/Lm0*PsiT; Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; Efmax=1.3889; Ven=1; delta=0; sd=1; ssd=0; Xd=1; Xdisp=0.2; Xq=1; n=-1; Ef=0.69:.002:.72; sq=1; ssq=0; nq=-1; d=0.01:0.01:pi; k=0; hold on for a=1:length(Ef) ef=Ef(a); while(abs(sd-ssd)>.0001 || (abs(sq-ssq)>.0001)) n=n+1; k=k+1; if n==0 %Condicional para no iterar en caso de no estar saturado sq=(sq+ssq); sd = (sd+ssd); end if n>0 %Condicional para el resto de las iteraciones

Page 8: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía sq=(sq+ssq)/2; sd=(sd+ssd)/2; end Xds = Xd/sd+((sd-1)/sd)*Xdisp; Xqs=Xq/sq+(sq-1)/sq*Xdisp; P(k)=((ef*sin(d(k)))/Xds)+(sin(2*(d(k)))*((1/Xqs)-(1/Xds))/2); Q(k)=ef/Xds*cos(d(k))-((cos(d(k))^2)/Xds+((sin(d(k))^2)/Xqs)); ie=(P(k)-1i*Q(k)); while(abs(sq-ssq)>.0001) Xqs=Xq/sq+((sq-1)/sq)*Xdisp; D=1+Xqs*1i*ie; Ee=1+Xdisp*1i*ie; Ed=abs(Ee)*sin(angle(D)-angle(Ee)); iflq=abs(Ed)/1.8; Psimq =abs(Ed); Mf=1/Lmsat; Mi=(1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsq=(Mf-Mi)/pi*(((Psimq-PsiT).*atan(tauT*(Psimq-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psimq-PsiT).^2)) ) + Psimq.*(Mf+Mi)/2; ssq=ifsq/iflq; sq=(sq+ssq)/2; end Eq=abs(Ee)*cos(angle(D)-angle(Ee)); ifl=abs(Eq)/1.8; tauT = fT / PsiT * Lm0/Lmsat; Psim =abs(Eq); ifs=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssd=ifs/ifl; end n=-1; k=0; Xds=Xd/sd+((sd-1)/sd)*Xdisp; Xqs=Xq/sq+(sq-1)/sq*Xdisp; ef=abs(D)+(Xds-Xqs)*abs(ie); ifm(a)=ef*ssd/1.8; P=ef/Xds*sin(d)+0.5*sin(2*d)*((1/Xqs)-(1/Xds)); plot(d*180/pi,P,'r'); sd=1; ssd=0; sq=1; ssq=0; end grid a=0:0.1:180; Pn=0.85*ones(length(a),1); plot(a,Pn,'--k')

Page 9: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

4.- El punto de operación a potencia nominal y factor de potencia unitario.

Comenzaremos con y procederemos a resolver el ejercicio a través

del método inverso hasta que converjan ambos grados de saturación, apoyándonos del

código de MATLAB 3 procederemos a realizar este proceso iterativo hasta conseguir los

valores finales en que convergen :

Código de MATLAB 3 hold off clear all clc Lm0=2; Lmsat=0.2; PsiT=.93; fT=1; iT = 1/Lm0*PsiT; Psimax =

4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; Efmax=1.3889; Ven=1; delta=0; sd=1; ssd=0; Xd=1; sq=1; ssq=0; Xdisp=0.2; Xq=1; n=-1; hold on while(abs(sq-ssq)>.0001) n = n+1; if n==0 %Condicional para no iterar en caso de no estar saturado sq=(sq+ssq); end if n>0 %Condicional para el resto de las iteraciones sq=(sq+ssq)/2; end Xqs= Xq/sq+((sq-1)/sq)*Xdisp; Ie = 0.85; D = 1 + j*Xqs*Ie; Ee = 1 + j*Xdisp*Ie; Ed = abs(Ee)*sin(angle(D)-angle(Ee)); ifl = abs(Ed)/1.8; Psim = abs(Ed); Mi = (1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifs = (Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+.5/tauT*(log(1+(tauT*PsiT)^2) - log(1+tauT^2*(Psim-

PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssq = ifs/ifl; scatter(n,sq,'b','d','filled') end D = 1 + j*Xqs*Ie; Ee = 1 + j*Xdisp*Ie; Eq = abs(Ee)*cos(angle(D)-angle(Ee)); ifld = abs(Eq)/1.8; Psim = abs(Eq); Mi = (1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi);

Page 10: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía ifsd = (Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+.5/tauT*(log(1+(tauT*PsiT)^2) - log(1+tauT^2*(Psim-

PsiT).^2)) ) + Psim.*(Mf+Mi)/2; sd = ifsd/ifld; Xds= Xd/sd+((sd-1)/sd)*Xdisp; Id=Ie*sin(angle(D)); Ef=abs(D)+(Xds-Xqs)*Id; ifop=Ef*sd/1.8; xlabel('Iteraciones') ylabel('Grado de saturación del eje directo (Sq)') grid

El grado de saturación del eje en cuadratura estabilizó en el valor de 1.018, con éste

obtenemos lugar del eje en cuadratura y por consiguiente el del eje directo:

0 1 2 3 4 5 6 71

1.005

1.01

1.015

1.02

1.025

Iteraciones

Gra

do d

e s

atu

ració

n d

el eje

directo

(S

q)

Page 11: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Con estos valores de impedancias podemos hallar el valor de la fuerza electromotriz

del campo, y con ésta el valor de la corriente de campo.

39.962

Page 12: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

5.- El punto de operación a potencia de 30 MW y corriente de campo máxima.

Resolviendo el método iterativo inverso y apoyándonos en el código de MATLAB 4

se obtuvo el siguiente punto de operación:

Sq Sd Xqs (p.u.) Xds (p.u.) Ef (p.u.) if (p.u.) Q (p.u.) Ie (p.u.) fp

1 2.2357 1 0.5578 1.2748 1.5834 0.4793 0.56545 57.957 0.53056 ind

Código de MATLAB 4 clear all clc Lm0=2;Lmsat=0.2;PsiT=.93;fT=1; iT = 1/Lm0*PsiT;Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax;tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; d=0:0.001:pi; sd=1; ssd=0; sq=1; ssq=0; n=-1; k=0; Xd=1; Xq=1; Pcalc=0; b=0; ef=1.5834*1.8; Xdisp=0.2; c=-1; while (abs(sq-ssq)>.0001 || abs(sd-ssd)>.0001) n=n+1; k=k+1; if n==0 %Condicional para no iterar en caso de no estar saturado sq=(sq+ssq); sd = (sd+ssd); end if n>0 %Condicional para el resto de las iteraciones sd=(sd+ssd)/2; end Xds=Xd/sd+((sd-1)/sd)*Xdisp; Xqs=Xq/sq+((sq-1)/sq)*Xdisp; ef=1.5834*1.8/sd; while (abs(Pcalc-0.3)>0.01) b=b+1;

1 1.5 2 2.5 3 3.5 41

1.5

2

2.5

Iteraciones

Gra

do d

e s

atu

racio

n d

el eje

en c

uadra

tura

(S

d)

Page 13: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía Pcalc=(ef*sin(d(b)))/Xds+0.5*sin(2*d(b))*((1/Xqs)-(1/Xds)); end d(b); Q=((ef*cos(d(b)))/Xds)-((cos(d(b))^2/Xds)+(sin(d(b))^2/Xqs)); ie=0.3-1i*Q; while (abs(sq-ssq)>.0001) c=c+1; if c==0 %Condicional para no iterar en caso de no estar saturado sd = (sd+ssd); end if c>0 %Condicional para el resto de las iteraciones sq=(sq+ssq)/2; sd=(sd+ssd)/2; end Xqs=Xq/sq+((sq-1)/sq)*Xdisp; D=1+(1i*Xqs*ie); Ee=1+(1i*Xdisp*ie); Ed=abs(Ee)*sin(abs(angle(D)-abs(angle(Ee)))); iflq=abs(Ed)/1.8; tauT = fT / PsiT * Lm0/Lmsat; Psim=abs(Ed); Mf=1/Lmsat; Mi=(1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsq=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssq=ifsq/iflq; if ssq<=1 ssq=1; end end c=-1; Xqs=Xq/sq+((sq-1)/sq)*Xdisp; D=1+(1i*Xqs*ie); Eq=abs(Ee)*cos(abs(angle(D))-abs(angle(Ee))); ifld=abs(Eq)/1.8; Psim=abs(Eq); Mi=(1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsd=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssd=ifsd/ifld; if ssd<=1 ssd=1; end b=0;Pcalc=0; hold on scatter(k,sd,'b','d','filled') end D=1+(1i*Xqs*ie); Id=abs(ie)*sin(abs(angle(D))+abs(angle(ie))); Ef=abs(D)+(Xds-Xqs)*Id; Ifop=Ef*sd/1.8; grid on

Page 14: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

6.- El punto de operación a potencia de -40 MW y corriente de campo nominal.

Utilizando el código de MATLAB 4 se halló el punto de operación con las

especificaciones pedidas en la pregunta actual. Resultó que el eje directo no satura, por lo

tanto Sq=1. El valor final en que estableció el grado de saturación del eje en cuadratura fue:

Sq Sd Xqs (p.u.) Xds (p.u.) Ef (p.u.) if (p.u.) Q (p.u.) Ie (p.u.) fp

1 2.0749 1 0.5856 0.8428 1 -0.2446 0.46886 148.55 0.85313 ind

Código de MATLAB 4 clear all clc Lm0=2;Lmsat=0.2;PsiT=.93;fT=1; iT = 1/Lm0*PsiT;Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax;tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; d=-pi/2:0.001:pi/2; sd=1; ssd=0; sq=1; ssq=0; n=-1; k=0; Xd=1; Xq=1; Pcalc=0; b=0; ef=1.5834*1.8; Xdisp=0.2; c=-1; while (abs(sq-ssq)>.0001 || abs(sd-ssd)>.0001) n=n+1; k=k+1; if n==0 %Condicional para no iterar en caso de no estar saturado sq=(sq+ssq); sd = (sd+ssd); end

1 1.5 2 2.5 3 3.5 4 4.5 5 5.5 61

1.2

1.4

1.6

1.8

2

2.2

2.4

Iteraciones

Gra

do d

e S

atu

ració

n d

el eje

en c

uadra

tura

(S

d)

Page 15: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía if n>0 %Condicional para el resto de las iteraciones sd=(sd+ssd)/2; end Xds=Xd/sd+((sd-1)/sd)*Xdisp; Xqs=Xq/sq+((sq-1)/sq)*Xdisp; ef=1*1.8/sd; while (abs(Pcalc+0.4)>0.01) b=b+1; Pcalc=(ef*sin(d(b)))/Xds+0.5*sin(2*d(b))*((1/Xqs)-(1/Xds)); end d(b); Q=((ef*cos(d(b)))/Xds)-((cos(d(b))^2/Xds)+(sin(d(b))^2/Xqs)); ie=-0.4-1i*Q; while (abs(sq-ssq)>.0001) c=c+1; if c==0 %Condicional para no iterar en caso de no estar saturado sd = (sd+ssd); end if c>0 %Condicional para el resto de las iteraciones sq=(sq+ssq)/2; sd=(sd+ssd)/2; end Xqs=Xq/sq+((sq-1)/sq)*Xdisp; D=1+(1i*Xqs*ie); Ee=1+(1i*Xdisp*ie); Ed=abs(Ee)*sin(abs(angle(D)-abs(angle(Ee)))); iflq=abs(Ed)/1.8; tauT = fT / PsiT * Lm0/Lmsat; Psim=abs(Ed); Mf=1/Lmsat; Mi=(1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsq=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssq=ifsq/iflq; if ssq<=1 ssq=1; end end c=-1; Xqs=Xq/sq+((sq-1)/sq)*Xdisp; D=1+(1i*Xqs*ie); Eq=abs(Ee)*cos(abs(angle(D))-abs(angle(Ee))); ifld=abs(Eq)/1.8; Psim=abs(Eq); Mi=(1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsd=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssd=ifsd/ifld; if ssd<=1 ssd=1;

Page 16: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía end b=0;Pcalc=0; hold on scatter(k,sd,'b','d','filled') end D=1+(1i*Xqs*ie); Id=abs(ie)*sin(abs(angle(D))+abs(angle(ie))); Ef=abs(D)+(Xds-Xqs)*Id; Ifop=Ef*sd/1.8; grid on xlabel('Iteraciones') ylabel('Grado de Saturación del eje en cuadratura (Sq)')

Page 17: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

7.- La característica de potencia activa en función del ángulo de carga.

Como se puede observar en la siguiente gráfica, la saturación en la máquina afecta

de manera directa la capacidad de generar potencia activa de la máquina, hasta el punto de

disminuir en un 15% la potencia máxima que es capaz de generar. Además se puede

observa un corrimiento de la gráfica, el cual produce que la potencia máxima generada sea

con un ángulo de carga (delta) cercano a 110°

Código de MATLAB 5

hold off clear all clc Lm0=2; Lmsat=0.2; PsiT=.93; fT=1; iT = 1/Lm0*PsiT; Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; s=2.0532; ss=0; Xd=1; Xdisp=0.2; Xq=1; n=-1; sq=1; ssq=0; nq=-1; d=0.001:0.01:pi; k=0; hold on; ifop=1.5834; for a=1:length(d) ef(a)=1.8*ifop/s; D=d(a); Xds=Xd/s+(s-1)/s*Xdisp; Xqs=Xq/sq+(sq-1)/sq*Xdisp; P(a)=((ef(a)*sin(D))/Xds)+(sin(2*(D))*((1/Xqs)-(1/Xds))/2); Q(a)=(ef(a)/Xds)*cos(D)-((cos(D)^2)/Xds+((sin(D)^2)/Xqs));

0 20 40 60 80 100 120 140 160 1800

0.5

1

1.5

2

2.5

3

Angulo de carga

P

Característica con Saturación

Característica sin Saturación

Page 18: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía ie(a)=(P(a)-1i*Q(a)); while((abs(sq-ssq)>.0001) && (abs(s-ss)>.0001)) Xqs=Xq/sq+((sq-1)/sq)*Xdisp; De=1+Xqs*1i*ie; Ee=1+Xdisp*1i*ie; Ed=abs(Ee)*sin(angle(De)-angle(Ee)); iflq=abs(Ed)/1.8; Psimq =abs(Ed); Mf=1/Lmsat; Mi=(1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsq=(Mf-Mi)/pi*(((Psimq-PsiT).*atan(tauT*(Psimq-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psimq-PsiT).^2)) ) + Psimq.*(Mf+Mi)/2; ssq=ifsq/iflq; if ssq<1 ssq=1; end Xqs=Xq/ssq+((ssq-1)/ssq)*Xdisp; De=1+Xqs*1i*ie; Eq=abs(Ee)*cos(angle(De)-angle(Ee)); ifl=abs(Eq)/1.8; tauT = fT / PsiT * Lm0/Lmsat; Psim =abs(Eq); ifs=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ss=ifs/ifl; if ss<1 ss=1; end sq=(sq+ssq)/2; s=(s+ss)/2; end Xds=Xd/s+(s-1)/s*Xdisp; Xqs=Xq/sq+(sq-1)/sq*Xdisp; P2(a)=ef(a)/Xds*sin(D)+0.5*sin(2*D)*((1/Xqs)-(1/Xds)); end % plot(d*180/pi,P,'r'); plot(d*180/pi,P2) grid on xlabel('Angulo de carga') ylabel('P')

Page 19: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

8.- El lugar geométrico de la corriente de armadura que no viola límites de operación.

En este punto procederemos a hallar las corrientes de armaduras necesarias para los

distintos factores de potencia tales que no se viole la corriente de campo máxima.

Código en MATLAB 6

Lm0=2;Lmsat=0.2;PsiT=.93;fT=1; iT = 1/Lm0*PsiT;Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax;tauT = fT / PsiT * Lm0/Lmsat; Mf = 1/Lmsat; fin=acos(.85); fi=fin:1*pi/180:pi/2; mod=1; Xdisp=.2; Xqs=1; sd=1;ssd=0;sq=1;ssq=0;n=-1;t=0;ifmax=2; Xq=1;Xd=1; for k=1:length(fi); while ifmax>=1.5802 while (abs(sq-ssq)>.0001) n=n+1; t=t+1; if n==0 %Condicional para no iterar en caso de no estar saturado sq=(sq+ssq); end if n>0 %Condicional para el resto de las iteraciones sq=(sq+ssq)/2; end Xqs=Xq/sq+((sq-1)/sq)*Xdisp; Ien=mod*(cos(fi(k))-1i*sin(fi(k))); Ee=1+1i*Ien*Xdisp; D=1+1i*Ien*Xqs; Ed=abs(Ee)*sin(angle(D)-angle(Ee)); ifl=abs(Ed)/1.8; tauT = fT / PsiT * Lm0/Lmsat; Psim=abs(Ed); Mf=1/Lmsat; Mi=(1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi);

0.2

0.4

0.6

0.8

1

30

210

60

240

90

270

120

300

150

330

180 0

Page 20: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía ifs=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ssq=ifs/ifl; if ssq<=1 ssq=1; end end Ien=mod*(cos(fi(k))-1i*sin(fi(k))); Ee=1+1i*Ien*Xdisp; D=1+1i*Xqs*Ien; Eq=abs(Ee)*cos(angle(D)-angle(Ee)); ifl=abs(Eq)/1.8; Psim=abs(Eq); Mi=(1/Lm0 - 1/Lmsat*(.5-atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifs=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; sd=ifs/ifl; Xqs=Xq/sq+((sq-1)/sq)*Xdisp; Xds=Xd/sd+((sd-1)/sd)*Xdisp; Df=1+1i*Xqs*Ien; Id=abs(Ien)*sin(angle(D)+fi(k)); Ef=abs(Df)+(Xds-Xqs)*Id; ifmax=Ef*sd/1.8; mod=mod-0.01; end end

Page 21: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

9.- Determine el triángulo de Potier de esta máquina.

En este punto, como no tenemos la característica de la corriente de cortocircuito en

función de la corriente de campo, supondremos por motivos didácticos que esta

característica es lineal y de pendiente igual a 1, a su vez supondremos una corriente de

armadura nominal para hallar la caída nominal en la reactancia de dispersión, la cual nos

dará la altura fija del triángulo de Potier (h=0.2). Con estas premisas junto a la

característica de vacío de la máquina se obtuvo el Triángulo de Potier y la característica de

saturación bajo carga:

Page 22: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

10.- Determine las curvas en V a tensión nominal para P=[0, 0.2, 0.4, 0.6, 0.8 y 1.0] pu.

Apoyándonos en el código de MATLAB 7, pudimos obtener las curvas en V de la

máquina sincrónica bajo estudio considerando la saturación de la misma (Gráfica morada),

donde se puede observar las corrientes de campo mínimas para obtener los distintos niveles

de potencia activa, además al compararla con la característica sin saturación (Gráfica azul)

podemos observar claramente como la saturación afecta considerablemente los distintos

puntos de operación de la máquina sincrónica

Page 23: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía

Código de MATLAB 7 clear all clc Lm0=2; Lmsat=0.2; PsiT=.93; fT=1; iT = 1/Lm0*PsiT; Psimax = 4*iT*Lmsat+PsiT; Psim = [0:0.002:1]*Psimax; tauT = fT / PsiT *

Lm0/Lmsat; Mf = 1/Lmsat; s=1; ss=0; Xd=1; Xdisp=0.2; Xq=1; n=-1; sq=1; ssq=0; nq=-1; d=0:0.001:pi; P=[1]; k=0; ifop=0.7:0.001:1.5834; ve=1; Pcalc=0; c=-1; for a=1:length(P) p=P(a); Xds=Xd/s+(s-1)/s*Xdisp; Xqs=Xq/sq+(sq-1)/sq*Xdisp; for b=1:length(ifop) if p==0 delta=0; end If=ifop(b); ef=If*1.8/s; if p>0 while (abs(Pcalc-p)>0.01) k=k+1; Pcalc=((ef/Xds)*sin(d(k)))+0.5*sin(2*d(k))*((1/Xqs)-

(1/Xds));

Page 24: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía end delta=d(k); end Q(b)=(ef/Xds)*cos(delta)-

((cos(delta)^2)/Xds+((sin(delta)^2)/Xqs)); ie(b)=(p-1i*Q(b)); while(abs(sq-ssq)>.0001 || (abs(s-ss)>.0001)) while(abs(sq-ssq)>.0001) c=c+1; if c==0 %Condicional para no iterar en caso de no estar

saturado sq = (sq+ss); end if c>0 %Condicional para el resto de las iteraciones sq=(sq+ssq)/2; end Xqs=Xq/sq+((sq-1)/sq)*Xdisp; De=1+Xqs*1i*ie; Ee=1+Xdisp*1i*ie; Ed=abs(Ee)*sin(angle(De)-angle(Ee)); iflq=abs(Ed)/1.8; if iflq==0 iflq=0.0001; end Psimq =abs(Ed); Mf=1/Lmsat; Mi=(1/Lm0 - 1/Lmsat*(.5-

atan(tauT*PsiT)/pi))/(.5+atan(tauT*PsiT)/pi); ifsq=(Mf-Mi)/pi*(((Psimq-PsiT).*atan(tauT*(Psimq-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psimq-PsiT).^2)) ) + Psimq.*(Mf+Mi)/2; ssq=ifsq/iflq; if ssq<1 ssq=1; end end c=-1; Xqs=Xq/sq+((sq-1)/sq)*Xdisp; De=1+Xqs*1i*ie; Eq=abs(Ee)*cos(angle(De)-angle(Ee)); ifl=abs(Eq)/1.8; tauT = fT / PsiT * Lm0/Lmsat; Psim =abs(Eq); ifs=(Mf-Mi)/pi*(((Psim-PsiT).*atan(tauT*(Psim-PsiT)) -

PsiT*atan(tauT*PsiT))+ .5/tauT*(log(1+(tauT*PsiT)^2) -

log(1+tauT^2*(Psim-PsiT).^2)) ) + Psim.*(Mf+Mi)/2; ss=ifs/ifl; if ss<1 ss=1; end s=(s+ss)/2; end Xds=Xd/s+(s-1)/s*Xdisp;

Page 25: Tarea N° 2. - prof.usb.veprof.usb.ve/jaller/calificaciones/Tarea2_CT4311_Macias.pdf · Universidad Simón Bolívar Departamento de Conversión y Transporte de Energía CT-4311 Conversión

Universidad Simón Bolívar

Departamento de Conversión y Transporte de Energía Xqs=Xq/sq+(sq-1)/sq*Xdisp; ef(a)=If*1.8/s; Q(b)=(ef(a)/Xds)*cos(delta)-

((cos(delta)^2)/Xds+((sin(delta)^2)/Xqs)); k=0; Pcalc=0; ie(b)=(p-1i*Q(b)); end end plot(ifop,abs(ie)) grid on xlabel('If') ylabel('Ie')