Se amplian los conceptos de circuitos pasivos analizados con transformadas e Laplace, aplicando a circuitos activos. Se obtienen los circuitos equivalentes o modelos matemáticos y se repiten los procedimientos anteriores.
El elemento activo más conocido es el amplificador operacional (op amp) que tienen ganancia «muy alta». El voltaje de salida v2 = – A v1, donde A es del orden de 105 o 106. Un factor importante es que la impedancia de entrada es muy alta del orden de 1012Ω y la impedancia de salida es muy baja (50-100Ω)
La configuración de la ganancia se establece con los resistores Ra y Rb y la forma de conectar las entradas y salidas.
Dada la alta impedancia del op amp, la corriente de retroalimentación I(s) fluye solo por los resistores. El voltaje de entrada es cero o muy pequeño dada la ganancia muy grande del op amp. Dadas estas simplificaciones, se aproxima con mucha precisión que:
1,1 Amplificador Operacional como multiplicador escalar
H(s)=−RRf
1,2 Amplificador Operacional como Integrador
Referencia: Oppenheim 11.52 p898
H(s)=(−RC1)s1
1,3 Amplificador Operacional como Sumador
Y(s)=−[R1RfX1(s)+R2RfX2(s)+R3RfX3(s)]
Observe que las ganancias del sumador son siempre negativas, hay una inversión de signo en cada señal de entrada.
Y(s)=K1X1(s)+K2X2(s)+K3X3(s)
Ejemplo 1. Implementación con Op-Amp
Referencia: Lathi 4.25 p401
Realizar la implementación de un sistema dado por la función de transferencia H(s)
H(s)=s2+4s+102s+5
El diagrama de bloques de la función de transferencia H(s) es,
Se agrupan algunos elementos como sumadores y sus factores de multiplicación. Para referencia se etiqueta cada punto como señal W(s) en cada punto donde el orden del exponente de ‘s’ es diferente.
Se considera la inversión de signo de la señal de entrada por la configuración del amplificador operacional y el factor K de cada rama a usar.
Se identifica el tipo de op amp a usar y se establecen los valores de los resistores en múltiplos de 10KΩ y los capacitores en el orden de 10 μF.
Ejemplos de solución a circuitos eléctricos que en el planteamiento se usan ecuaciones diferenciales ordinarias y se desarrollan usando la Transformada de Laplace. En el ejercicio 3 se desarrolla una malla que genera un sistema de ecuaciones EDO y se resuelve con transformadas de Laplace.
Ejercicio 1. Circuito RLC, fuente DC e interruptor con transformada de Laplace
Referencia: Lathi Ej.4.13 p365
En el circuito de la figura, el interruptor se encuentra cerrado mucho tiempo antes de t=0. Cuando se abre en un instante, encuentre la corriente del inductor y(t) para t≥0.
El interruptor se encuentra cerrado por mucho tiempo, la corriente por el inductor es 2 A y el voltaje del capacitor es 10 V, pues el inductor en DC opera como un conductor sin resistencia y el capacitor se encuentra completamente cargado.
Cuando se abre el interruptor, el circuito se equivalente es el mostrado en la derecha, La corriente inicial del inductor es y(0–)=2 y el voltaje inicial del capacitor Vc(0–)=10.
El voltaje en la entrada es 10 V, empezando en t=0 y puede ser representado como 10 μ(t) para simplificar la representación de la fuente DC y expresar que antes de t=0 se aplica las condiciones iniciales dadas.mientras se den las condiciones iniciales en t=0 solo será necesario saber que la corriente en t≥0 para determinar la respuesta en t≥0 .
La ecuación del circuito en luego de abrir el interruptor es:
Ejemplo 2. Circuito RLC con transformada de Laplace
Referencia: Lathi Ejemplo 4.17 p375
Realice el análisis de la corriente i(t) en el circuito mostrado en la figura, suponiendo que todas las condicioes iniciales son cero.
El primer paso es representar el circuito en el dominio de la frecuencia (s), mostrado a la derecha de la figura. Se representan los voltajes y corrientes usando la transformada de Laplace.
El voltaje 10μ(t) se representa como 10/s y la corriente i(t) como I(s). Todos los elementos del circuito se muestran con su respectiva impedancia. El inductor de 1 henrio se representa por s, el capacitor de 1/3 faradio se representa por 2/s. El resistor de 3 ohms es solo 3.
El voltaje en el circuito V(s) en cualquier elemento es I(s) por su impedancia Z(s).
La impedancia del circuito para el lazo es:
Z(s)=s+3+s2=ss2+3s+2
el voltaje de entrada es V(s)=10/s, por lo que la fórmula básica de corriente es:
Ejemplo 3. Circuito RC y RL, fuente DC e interruptor con transformada de Laplace
Referencia: Lathi ejemplo 4.18 p378
El interruptor está en posición cerrada por un largo tiempo antes de que sea abierto en t=0. Encuentre las corrientes y1(t) y y2(t) para t>0
Al revisar el circuito se observa que cuando el interruptor esta cerrado y se han alcanzado las condiciones de estado estable, el voltaje del capacitor Vc=16 V y la corriente del inductor y2= 4 .
Las fuentes de 20V y 4V se pueden también ordenar en forma equivalente
Cuando se abre el interruptor en t=0, las condiciones iniciales son Vc(0–)=16 y y2(0–)= 4 como se muestra en la versión con transformada de Laplace en la figura anterior. La ecuación diferencial del circuito para suma de voltajes en la malla en el lazo 1 y la tabla de propiedades de la transformada de Laplace para la integración aplicada al capacitor,
El voltaje inicial del capacitor se puede representar por un voltaje en serie 16/s y la corriente inicial del inductor de 4 A que representa una fuente de valor VL = Ly2(0–) = 1/2(4) = 2. La ecuación diferencial del lazo 2 y usando la transformada de Laplace para la derivada con condiciones iniciales para el inductor,
Los ejercicios incluyen constantes tipo símbolo, que por simplicidad de la ecuación en la resolución se usan solo las operaciones básicas para la transformada Inversa de Laplace.
Ejemplo 1. Y(s) dado H(s) y X(s) con una constante α
Referencia: Oppenheim ejemplo 9.37 p718 pdf749, semejante al ejercicio sistema-modelo LTI Lathi ejemplo 1.16 p111. Oppenheim problema 2.61c p164 pdf195
Suponga un sistema LTI causal que se describe mediante la ecuación diferencial:
Se crea la expresión H(s) y X(s) en forma simbólica usando las variables, s y a.
Dado que la expresión es tipo polinomio, y que no existen exponenciales, el desarrollo del algoritmo puede realizarse con operaciones básicas realizadas en las secciones anteriores. Añadiendo la operación para encontrar la salida Y(s)
donde Y(s) es la transformada unilateral de Laplace de y(t)
el último término, representa la respuesta del sistema LTI causal y la condición de reposo inicial, conocido también como la respuesta en estado cero.
Una interpretación análoga se aplica a los primeros dos términos del miembro derecho de la ecuación, que representa la respuesta del sistema cuando la entrada es cero (α=0) , conocida también como respuesta a entrada cero.
Observe que la respuesta a entrada cero es una función lineal de los valores de las condiciones iniciales. Demostrando que la respuesta es la superposición del estado cero y respuesta a entrada cero.
si α=2, β = 3 y γ=-5, al realizar la expansión en fracciones parciales encontramos que:
Y(s)=s1−s+11+s+23
con la aplicación de las tablas de transformadas se obtiene:
y(t)=[1−e−t+3e−2t]μ(t), para t>0
Instrucciones con Sympy-Python
En este caso, se ingresan las propiedades de diferenciación para Laplace, aunque pueden crearse a partir de una fórmula mas general en la siguiente sección.
La ecuación se define como una expresión de (lado_izquierdo)=(lado_derecho) contruida a partir la ecuación del problema incluyendo las condiciones iniciales diferentes de cero.
Luego se busca la expresión solución para Ys, por simplicidad se la cambia a fracciones parciales, asi es más sencillo usar la tabla de transformadas de Laplace.
Para una sustituir los valores de a,b y g, se usa un diccionario {} con las parejas de {variable: valor}.
ecuacion:
2 a
Y*s + 3*Y*s + 2*Y - b*s - 3*b - g = -
s
Y(s):
2
a + b*s + 3*b*s + g*s
----------------------
/ 2 \
s*\s + 3*s + 2/
Y(s) en Fracciones parciales
a a - 2*b - 2*g a - 2*b - g
--- + ------------- - -----------
2*s 2*(s + 2) s + 1
sustituye valores de a, b y g
Ys:
3 1 1
----- - ----- + -
s + 2 s + 1 s
y(t) con constantes simbolo
/ 2*t t\ -2*t
\a*e + a - 2*b - 2*g + 2*(-a + 2*b + g)*e /*e *Heaviside(t)
-----------------------------------------------------------------
2
sustituye valores de a, b y g
y(t):
/ 2*t t \ -2*t
\e - e + 3/*e *Heaviside(t)
y(t) expande expresion
-t -2*t
Heaviside(t) - e *Heaviside(t) + 3*e *Heaviside(t)
>>>
Compare los resultados con lo expuesto en la parte teórica.
# Y(s) con condiciones iniciales como constantes# Oppenheim ejemplo 9.38 p719import numpy as np
import sympy as sym
# INGRESO
t = sym.Symbol('t', real=True)
s = sym.Symbol('s')
a,b,g = sym.symbols('a b g')
Y = sym.Symbol('Y')
# propiedad de diferenciacion
d2Y = (s**2)*Y - b*s - g
dY = s*Y - b
Xs = a/s # señal de entrada# PROCEDIMIENTO# ecuacion del sistema
Lizq = d2Y + 3*dY + 2*Y # Lado izquierdo
Lder = Xs
ecuacion = sym.Eq(Lizq, Lder)
# despeja la ecuacion para Y(s)
Ys = sym.solve(ecuacion,Y)[0]
Yspc = sym.apart(Ys,s) # incluye constante simbolo
Ysp = Yspc.subs({a:2,b:3,g:-5})
# y(t) con inversa transformada Laplace
ytc = sym.inverse_laplace_transform(Yspc,s,t)
yt = sym.inverse_laplace_transform(Ysp,s,t)
# SALIDAprint('ecuacion:')
sym.pprint(ecuacion)
print('\n Y(s): ')
sym.pprint(Ys)
print('\n Y(s) en Fracciones parciales')
sym.pprint(Yspc)
print('\nsustituye valores de a, b y g')
print(' Ys:')
sym.pprint(Ysp)
print('\n y(t) con constantes simbolo ')
sym.pprint(ytc)
print('\nsustituye valores de a, b y g')
print('\n y(t): ')
sym.pprint(yt)
print('\n y(t) expande expresion')
sym.pprint(sym.expand(yt))
Se presenta un ejercicio que integra los temas tratados en la Unidad 4 para el análisis del Sistema con Transformada de Laplace. El algoritmo se muestra como una referencia para adaptar en los ejercicios integradores usando las funciones resumidas en telg1001.py .
Referencia: Lathi Ejemplo 4.12 p361
Resuelva la ecuación diferencial de segundo orden
(D2+5D+6)y(t)=(D+1)x(t)
para las condiciones iniciales y(0-)=2 y y'(0-)=1, cuando se presenta la entrada
x(t)=e−4tμ(t)
Desarrollo
La ecuación diferencial del sistema, interpretando los operadores D por derivadas se muestra a continuación.
δt2δ2y(t)+5δtδy(t)+6y(t)=δt2δ2x(t)+x(t)
La ecuación se representa en el dominio s de las transformadas de Laplace para usarlas en el algoritmo:
con lo que se escribe la función de transferencia H(s) que se ingresa para el análisis del algoritmo.
H(s)=Y(s)X(s)=s2+5s+6s+1
Para el caso de x(t) se usa la transformada de Laplace con el algoritmo para obtener X(s), mostrando que se puede realizar directamente para términos simples. Si la expresión es mas elaborada se usaría las funciones realizadas para f(t) con Sympy-Python
Con el algoritmo se obtienen los siguientes resultados, con lo que se puede analizar la estabilidad del sistema usando la gráfica de polos. No se observan polos en el lado derecho del plano RHP, por lo que el sistema es asintóticamente estable, BIBO-estable.
la gráfica de H(s) junto a los polos:
la gráfica de h(t) como respuesta al impulso,
la gráfica de la señal de entrada x(t) y la respuesta y(t) es:
Las expresiones para Y(s), H(s), análisis de estabilidad que entrega el algoritmo son:
Con los datos de H(s) y H(s) con fracciones parciales, se presentan dos formas de realizar el diagrama de bloques. El primer diagrama corresponde a H(s), mientras qu el segundo corresponde a la interpretación de las fracciones parciales.
Instrucciones Python
Usando los bloques desarrollados en la Unidad 4 Sistemas LTI – Laplace y las funciones resumidas como telg1001.py que pueden ser usados en cada pregunta, se agrupan y tiene como resultado:
# Y(s) Respuesta total con entada cero y estado cero# Qs Y(s) = Ps X(s) ; H(s)=Ps/Qs# https://blog.espol.edu.ec/telg1001/import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
# INGRESO
s = sym.Symbol('s')
t = sym.Symbol('t', real=True)
d = sym.DiracDelta(t)
u = sym.Heaviside(t)
# H(s) respuesta impulso
Ps = s+1
Qs = s**2 + 5*s + 6
Hs = Ps/Qs
#Hs = 1+0*s cuando es constante# X(s) Señal de entrada
xt = sym.exp(-4*t)*u
# condiciones iniciales, [y'(0),y(0)] orden descendente
t0 = 0
cond_inicio = [1, 2] # estado cero no se usan# Grafica, intervalo tiempo [t_a,t_b]
t_a = -1 ; t_b = 10
muestras = 101 # 51 resolucion grafica# PROCEDIMIENTO
Hs = fcnm.apart_s(Hs) # fracciones parciales
Hs_fc = fcnm.factor_exp(Hs) # en factores
Hs_Qs2 = fcnm.Q_cuad_s_parametros(Hs_fc)
polosceros = fcnm.busca_polosceros(Hs)
Q_polos = polosceros['Q_polos']
P_ceros = polosceros['P_ceros']
estable = fcnm.estabilidad_asintotica_s(Q_polos)
# H(t) respuesta al impulso
ht = 0*s
term_suma = sym.Add.make_args(Hs)
for term_k in term_suma:
ht_k = sym.inverse_laplace_transform(term_k,s,t)
# simplifica log(exp()) ej: e**(-2s)/(s**2)if ht_k.has(sym.log):
ht_k = sym.simplify(ht_k,inverse=True)
ht = ht + ht_k
lista_escalon = ht.atoms(sym.Heaviside)
ht = sym.expand(ht,t) # terminos suma
ht = sym.collect(ht,lista_escalon)
# PROCEDIMIENTO Respuesta ZIR, ZSR
Xs = fcnm.laplace_transform_suma(xt)
# ZIR_s respuesta entrada cero de s
sol_ZIR = fcnm.respuesta_ZIR_s(Hs,cond_inicio)
ZIR = sol_ZIR['ZIR']
yt_ZIR = sol_ZIR['yt_ZIR']
# ZSR respuesta estado cero, Y(s) a entrada X(s)
sol_ZSR = fcnm.respuesta_ZSR_s(Hs,Xs)
ZSR = sol_ZSR['ZSR']
yt_ZSR = sol_ZSR['yt_ZSR']
# Respuesta total Y(s) y y(t)
Ys = ZIR + ZSR
Ys = fcnm.apart_s(Ys)
yt = yt_ZIR + yt_ZSR
lista_escalon = yt.atoms(sym.Heaviside)
yt = sym.collect(yt,lista_escalon)
# SALIDAprint(' H(s) = P(s)/Q(s):')
sym.pprint(Hs)
print(' H(s) en factores:')
sym.pprint(Hs_fc)
iflen(Hs_Qs2)>0:
print('\nH(s) parámetros cuadraticos:')
fcnm.print_resultado_dict(Hs_Qs2)
print('\n h(t) :')
sym.pprint(ht)
print('\npolosceros:')
fcnm.print_resultado_dict(polosceros)
print('\nEstabilidad de H(s):')
for k in estable:
print('',k,':',estable[k])
print('\n X(s): ')
sym.pprint(Xs)
print('\nRespuesta entrada cero ZIR H(s) y condiciones iniciales')
ifnot(sol_ZIR == sym.nan): # existe resultado
fcnm.print_resultado_dict(sol_ZIR)
else:
print(' insuficientes condiciones iniciales')
print(' revisar los valores de cond_inicio[]')
print('\n ZSR respuesta estado cero:')
fcnm.print_resultado_dict(sol_ZSR)
print('\n Y(s)_total = ZIR + ZSR:')
sym.pprint(Ys)
print('\n y(t)_total = ZIR + ZSR:')
sym.pprint(yt)
# Graficas polos, H(s), con polos h(t) --------
figura_s = fcnm.graficar_Fs(Hs,Q_polos,P_ceros,f_nombre='H',solopolos=True)
figura_Hs = fcnm.graficar_Fs(Hs,Q_polos,P_ceros,muestras,f_nombre='H')
figura_ht = fcnm.graficar_ft(ht,t_a,t_b,muestras,f_nombre='h')
# GRAFICAS y(t),x(t),h(t) ---------------------
figura_ft = fcnm.graficar_xh_y(xt,ht,yt,t_a,t_b,muestras)
plt.show()
1. Respuesta a entrada cero y condiciones iniciales
Si el sistema bajo observación se no tiene señal de entrada X(s) = 0, la parte de estado cero no aporta componente de salida.
Respuesta total
=
respuesta a entrada cero ZIR
+
respuesta a estado cero ZSR
Sin embargo, si existen valores iniciales en los componentes del sistema (circuito con capacitores), se producirá una respuesta a entrada cero que básicamente responde a las condiciones iniciales al sistema.
Referencia: Lathi Ejemplo 4.12 p361
Un sistema LTI tiene la ecuación diferencial de segundo orden expresada con operadores D:
(D2+5D+6)y(t)=(D+1)x(t)
El sistema tiene las condiciones iniciales y(0-)=2 y y'(0-)=1.
1.1 Desarrollo Analítico para Respuesta a entrada cero y condiciones iniciales
La ecuación se escribe usando la transformada de Laplace:
(s2+5s+6)Y(s)=(s+1)X(s)
El diagrama básico de la expresión LTI, sin simplificar, se muestra en la figura. La sección de color naranja representa Q(s) y la sección en verde representa P(s)
Para el análisis de entrada cero X(s)=0, la ecuación se simplifica a
(s2+5s+6)Y(s)=0
Las condiciones iniciales se aplican usando la tabla de propiedades de Laplace, donde se muestra que los valores iniciales se añaden a los componentes de primera y segunda derivada:
Se inicia con los resultados de algoritmos anteriores, donde se ingresan las expresiones de Q(s) y P(s) correspondientes a H(s). Las condiciones iniciales se ingresan en un arreglo en orden de derivada descendente, por ejemplo para el ejercicio [y'(0), y(0)] = [1., 2.]
Par generar las expresiones de condicion inicial, se genera el término de condición cero o term_cero, compuesto de los términos de orden tomados del vector inicial.
# Y_ZIR(s) Respuesta entrada cero con Transformada de Laplace# Qs Y(s) = Ps X(s) ; H(s)=Ps/Qs# Ejemplo Lathi 4.12 p361# https://blog.espol.edu.ec/telg1001/import sympy as sym
import telg1001 as fcnm
# INGRESO
s = sym.Symbol('s')
t = sym.Symbol('t', real=True)
# H(s) respuesta impulso
Ps = s + 1
Qs = s**2 + 5*s + 6
Hs = Ps/Qs
#Hs = 6*sym.exp(-2*s)/s**2# X(s) Señal de entrada
Xs = 0*s
# condiciones iniciales, [y'(0),y(0)] orden descendente
t0 = 0
cond_inicio = [1, 2]
# PROCEDIMIENTOdefrespuesta_ZIR_s(Hs,cond_inicio):
'''respuesta a entrada cero ZIR con H(s) y condiciones de inicio
Si ZIR es sym.nan existen insuficientes condiciones
'''
polosceros = fcnm.busca_polosceros(Hs)
Q_polos = polosceros['Q_polos']
# coeficientes y grado de Q
Q = 1+0*s
for raiz in Q_polos:
veces = Q_polos[raiz]
Q = Q*((s-raiz)**veces)
Q = sym.poly(Q,s)
Q_coef = Q.all_coeffs() # coeficientes Q
Q_grado = Q.degree() # grado polinomio Q
term_cero = 0*s
iflen(cond_inicio) == Q_grado:
for j inrange(0,Q_grado,1):
term_orden = 0
for i inrange(j,Q_grado,1):
term_orden = term_orden + cond_inicio[i]*(s**(i-j))
term_cero = term_cero + term_orden*Q_coef[j]
ZIR = term_cero/Q
ZIR = fcnm.apart_s(ZIR)
ZIR_Qs2 = fcnm.Q_cuad_s_parametros(ZIR)
else:
ZIR = sym.nan # insuficientes condiciones inicialesifnot(ZIR==sym.nan):
# entrada_cero en t
yt_ZIR = sym.inverse_laplace_transform(ZIR,s,t)
# simplifica log(exp()) ej: e**(-2t)/(s**2)if yt_ZIR.has(sym.log):
yt_ZIR = sym.simplify(yt_ZIR,inverse=True)
lista_escalon = yt_ZIR.atoms(sym.Heaviside)
yt_ZIR = sym.expand(yt_ZIR,t) # terminos suma
yt_ZIR = sym.collect(yt_ZIR,lista_escalon)
sol_ZIR ={'term_cero' : term_cero}
iflen(ZIR_Qs2)>0: # añade si Q es cuadratico
sol_ZIR['ZIR_Qs2'] = ZIR_Qs2
sol_ZIR = sol_ZIR | {'ZIR' : ZIR,
'yt_ZIR': yt_ZIR}
else:
sol_ZIR = sym.nan
return(sol_ZIR)
# respuesta a entrada_cero ZIR en s
sol_ZIR = respuesta_ZIR_s(Hs,cond_inicio)
# SALIDA - polinomioprint('\nRespuesta entrada cero ZIR H(s) y condiciones iniciales')
ifnot(sol_ZIR == sym.nan): # existe resultadofor k in sol_ZIR:
print('\n',k,':')
sym.pprint(sol_ZIR[k])
else:
print(' insuficientes condiciones iniciales')
print(' revisar los valores de cond_inicio[]')
El siguiente caso sería realizar el mismo sistema LTI CT, con respuesta Y(s) ante una extrada X(s), conocido también como estado cero.
Para mantenerlo simple, se consideraría condiciones iniciales en cero.
Finalmente se sumarían la entrada cero y estado cero para obtener la respuesta total.
2. Respuesta a estado cero con x(t) de Sistema con ecuación diferencial de 2do orden
El ejercicio continua el desarrollo del concepto sobre la respuesta total del sistema tiene dos componentes:
Respuesta total
=
respuesta a entrada cero
+
respuesta a estado cero
ZSR = h(t) ⊗ x(t)
ZSR(s) = H(S) X(s)
Referencia: Lathi Ejemplo 4.12 p361
Un sistema LTI tiene la ecuación diferencial de segundo orden expresada con operadores D:
(D2+5D+6)y(t)=(D+1)x(t)
(viene de la sección anterior…) si se somete a la entrada,
x(t)=e−4tμ(t)
2.1. Desarrollo Analítico para estado cero
Ahora, todos los componentes tienen condición inicial cero y una señal de entrada x(t). El diagrama de bloques del sistema se obtuvo en la sección anterior de «entrada cero»:
En este caso se considera que ningún componente del sistema tiene valor inicial [0., 0.], interpretado por ejemplo como que todos los capacitores estan descargados, la salida corresponde solo a la parte conformada por x(t) y la respuesta al impulso h(t).
para la entrada x(t),
x(t)=e−4tμ(t)⟺X(s)=s+41
Usando transformadas de Laplace, la ecuación del sistema se convierte en:
# Y_ZSR(s) Respuesta estado cero con Transformada Laplace# Qs Y(s) = Ps X(s) ; H(s)=Ps/Qs# Ejemplo Lathi 4.12 p361# https://blog.espol.edu.ec/telg1001/import sympy as sym
import telg1001 as fcnm
# INGRESO
s = sym.Symbol('s')
t = sym.Symbol('t', real=True)
# H(s) respuesta impulso
Ps = s + 1
Qs = s**2 + 5*s + 6
Hs = Ps/Qs
# X(s) Señal de entrada
Xs = 1/(s+4)
# PROCEDIMIENTOdefrespuesta_ZSR_s(Hs,Xs):
'''respuesta a estado cero ZSR con H(s) y X(s)
'''
ZSR = Hs*Xs
ZSR = fcnm.apart_s(ZSR)
ZSR_Qs2 = fcnm.Q_cuad_s_parametros(ZSR)
# ZSR respuesta estado cero(t)
yt_ZSR = sym.inverse_laplace_transform(ZSR,s,t)
# simplifica log(exp()) ej: e**(-2t)/(s**2)if yt_ZSR.has(sym.log):
yt_ZSR = sym.simplify(yt_ZSR,inverse=True)
lista_escalon = yt_ZSR.atoms(sym.Heaviside)
yt_ZSR = sym.expand(yt_ZSR,t) # terminos suma
yt_ZSR = sym.collect(yt_ZSR,lista_escalon)
sol_ZSR = {'ZSR' : ZSR}
iflen(ZSR_Qs2)>0:
sol_ZSR['ZSR_Qs2'] = ZSR_Qs2
sol_ZSR = sol_ZSR | {'yt_ZSR' : yt_ZSR}
return(sol_ZSR)
sol_ZSR = respuesta_ZSR_s(Hs,Xs)
# SALIDA - polinomioprint('\n ZSR respuesta estado cero:')
for k in sol_ZSR:
print('\n',k,':')
sym.pprint(sol_ZSR[k])
El siguiente bloque consiste en hacer el ejercicio completo del sistema LTIC, considerando condiciones iniciales y una entrada x(t), que da la respuesta total del sistema.
En un sistema acotado (bounded) estable, sometido a una señal de entrada en la que se conoce que su salida siempre es la misma, se lo considera un sistema BIBO-estable. (Bounded input, Boundes output)
La función de transferencia H(s) de un sistema es una descripción externa con la forma:
H(s)=Q(s)P(s)
El polinomio Q(s) representa una descripción interna, en la que si todos los polos de H(s) se encuentran en la parte izquierda del plano (left half plane, LHP), todos los términos h(t) son exponenciales decrecientes y h(t) es absolutamente integrable. En consecuencia el sistema es BIBO-estable.
Se asume que H(s) es una función en la que M ≤ N, siendo M el grado de P(s) y N e grado de Q(s), caso contrario el sistema es es BIBO-inestable.
1. Un sistema LTI CT es asintóticamente estable si y solo si todos los polos de su función de transferencia H(s) se encuentran en el lado izquierdo del plano, pudiendo ser simples o repetidos.
2. Un sistema LTI CT es inestable si y solo si se dan una o ambas de las condiciones:
2.1 Existe al menos un polo de H(s) en el lado derecho del plano RHP.
2.2 Existen polos repetidos de H(s) en el eje imaginario.
3. Un sistema LTI CT es marginalmente estable si y solo si, no hay polos de H(s) en el lado derecho del plano RHP y algunos polos NO repetidos sobre en el eje imaginario (parte Real = 0).
El grado del polinomio Q es N=3 y es mayor que el grado del polinomio P es M=1.
La expresión de entrada se escribe en el algoritmo como:
Hs = ((s-3)/(s+1))*(1/(s**2+4*s+8))
Los polos se obtienen como las raíces del denominador Q, necesario para plantear las fracciones parciales y simplificar la expresión como la suma de componentes más simples.
polos: {-1: 1, -2 - 2*I: 1, -2 + 2*I: 1}
Las fracciones parciales de H(s) usando el algoritmo de la sección anterior se obtiene
En la gráfica se observa que todos los polos se encuentran en el lado izquierdo del plano LHP. Por lo que se considera que existen asintotas en la función correspondiente a tiempo, y el sistema es «asintóticamente estable«. También se considera BIBO estable por no tener los polos en el lado derecho RHP=0
la gráfica de h(t) muestra lo indicado por el análisis de polos
Instrucciones en Python Ejemplo 1
El algoritmo es una continuación del manejo de Polos y Ceros realizado en la sección anterior, donde se añade la revisión de las raíces de Q cuando son de tipo imaginarias.
Por simplicidad para el análisis se usan solo las partes esenciales, asi el enfoque es sobre la posición de los polos y su relación con h(t).
En el proceso, se agrupan las partes reales en Q_real y las imaginaria en Q_imag para realizar el conteo de lo repetido.
# Polos y ceros de funcion de transferencia H(s)# Estabilidad del sistema, H(s) respuesta al impulso# Ps es numerador, Qs es denominadorimport numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
# INGRESO
s = sym.Symbol('s')
t = sym.Symbol('t',real=True)
Hs = ((s-3)/(s+1))*(1/(s**2+4*s+8))
# Grafica, intervalo tiempo [t_a,t_b]
t_a = -1 ; t_b = 10
muestras = 101 # 51 resolucion grafica# PROCEDIMIENTO
polosceros = fcnm.busca_polosceros(Hs)
P_ceros = polosceros['P_ceros']
Q_polos = polosceros['Q_polos']
# h(t) desde H_fp(s) con inversa Laplace
ht = sym.inverse_laplace_transform(Hs,s,t)
# SALIDA# Analiza estabilidad asintotica
casi_cero = 1e-12
# Separa parte real e imaginaria de raices
cuenta_real = 0; cuenta_imag = 0
unicos = 0 ; repetidos = 0 ; enRHP = 0
for raiz in Q_polos:
[r_real,r_imag] = raiz.as_real_imag()
ifabs(r_real)>casi_cero andabs(r_imag)<casi_cero :
cuenta_real = cuenta_real+1
# para estabilidad asintotica
conteo = Q_polos[raiz]
if conteo==1 and r_real==0 andabs(r_imag)>0:
unicos = unicos + 1
if conteo>1 and r_real==0 andabs(r_imag)>0:
repetidos = repetidos + 1
if r_real>0:
enRHP = enRHP + 1
cuenta_imag = len(Q_polos)-cuenta_real
# Revisa lado derecho del plano RHP
asintota = ""if enRHP==0:
asintota = 'estable'if enRHP>0 or repetidos>0:
asintota = 'inestable'if enRHP==0 and unicos>0:
asintota = 'marginalmente estable'# SALIDAprint('busca_PolosCeros:')
for k in polosceros.keys():
eq_lista = ['H','H_fp']
if k in eq_lista:
print('',k,':')
sym.pprint(polosceros[k])
else:
print('\n',k,':',polosceros[k])
print(' polos reales: ',cuenta_real)
print(' polos complejos: ',cuenta_imag)
print(' sobre lado derecho RHP:',enRHP)
print(' sobre Eje Imaginario, repetidos: ',
repetidos,' unicos:', unicos)
print('\n asintoticamente: ', asintota)
print('\n h(t): ')
sym.pprint(ht)
# GRAFICA -----------
muestras_H = 161
fig_Hs = fcnm.graficar_Fs(Hs,Q_polos,P_ceros,
muestras=muestras_H,
f_nombre='H',
solopolos=True)
figura_h = fcnm.graficar_ft(ht,t_a,t_b,f_nombre='h')
plt.show()
El grado del polinomio Q es N=3 y es mayor que el grado del polinomio P es M=1.
Los polos se obtienen como las raíces del denominador Q, necesario para plantear las fracciones parciales y simplificar la expresión en componentes más simples.
La expresión para el algoritmo se escribe como:
Hs = ((s+2)/(s-1))*(1/(s**2+4*s+8))
y el resultado de polos obtenido es:
polos: {1: 1, -2 - 2*I: 1, -2 + 2*I: 1}
Las fracciones parciales de H(s) usando el algoritmo de la sección anterior se obtiene
En la gráfica se observa que un polo se encuentra en el lado derecho del plano RHP. Un polo en RHP genera un componente creciente en el tiempo, en consecuencia, el sistema es «asintóticamente inestable». También el sistema es BIBO inestable.
la gráfica de h(t) complementa lo interpretado con la posición de los polos en el plano s
Ejemplo 3 H(s) con polos complejos sobre eje imaginario únicos
Referencia: Lathi 2.14c p201
Revise la estabilidad BIBO de los siguientes sistemas LTIC descritos por las ecuaciones:
(D+2)(D2+4)y(t)=(D2+D+1)x(t)
Se expresa la ecuación de operadores diferenciales D por la transformada de Laplace.
(s+2)(s2+4)Y(s)=(s2+s+1)X(s)
Ps = s**2+s+1
Qs = (s+2)*(s**2+4)
Hs = Ps/Qs
El grado del polinomio Q es N=3 y el grado del polinomo P es M=2, los polos del denominador muestran raíces complejas, por lo que se tendrá componentes cuadráticos para la transformada inversa de Laplace.
Q_polos : {-2: 1, -2*I: 1, 2*I: 1}
Usando fracciones parciales con el algoritmo se tiene que,
Ejemplo 4. H(s) con polos repetidos sobre eje imaginario
Referencia: Lathi 2.14d p201
Revise la estabilidad BIBO de los siguientes sistemas LTIC descritos por las ecuaciones:
(D+1)(D2+4)2y(t)=(D2+2D+8)x(t)
Se expresa la ecuación de operadores diferenciales D por la transformada de Laplace.
(s+1)(s2+4)2Y(s)=(s2+2s+8)X(s)
Las raíces del denominador Q muestran que tienen valores complejos, mostrando que se usaría un componente cuadratico para la transformada inversa de Laplace.
polos: {-1: 1, -2*I: 2, 2*I: 2}
al presentar dos raíces repetidas sobre el eje imaginario, se considera asintóticamente inestable y también BIBO inestable.
Revise el resultado de la inversa de Laplace usando las tablas.
Todos los sistemas de procesamiento de señales deben ser asintóticamente estables.
Los sistemas inestables no son requeridos debido a que con condiciones iniciales, intencionales o no, llevan a respuestas fuera de rango. Estas respuestas inestables destruyen el sistema o lo llevan a condiciones de saturación que cambian la naturaleza del sistema debido al tipo de crecimiento exponencial.
Sistemas marginalmente estables, aunque BIBO inestables, tienen aplicación como osciladores, que es un sistema que genera su propia señal sin la aplicación de entradas externas.
La salida de un oscilador es una respuesta de entrada cero. Si la respuesta es una sinusiode con frecuencia ω0 el sistema es marginalmente estable en ±jω0. Para diseñar un oscilador de frecuencia se debe usar un sistema descrito por la ecuación diferencial:
(D2+ω02)y(t)=x(t)
Sin embargo, osciladores en la práctica se construyen con sistemas no lineales.
Referencia: Lathi ejemplo 4 p330, Hsu literal C. p112
La función de transferencia H(s) se escribe como P(s)/Q(s) que es una expresión racional de s. Se considera que en la expresión el grado m del polinomio del numerador P(s) es mayor que el grado n del polinomio del denominador Q(s). Las raíces del polinomio P(s) del numerador se conocen como ceros, pues convierten la expresión en cero, mientras que para el denominador las raices se denominan polos al convertir en infinito la expresión (Hsu).
La expresión H(s) tiene los componentes de P(s) y Q(s) como se indican,
H(s)=s2+4s+32s+4
Ingreso de H(s) como P(s) y Q(s)
Para la expresión de la función de transferencia, con Sympy se crean los polinomios del numerador Ps y denominador Qs.
# INGRESO
s = sym.Symbol('s')
Ps = 2*s + 4 # 1+0*s cuando es constante
Qs = s**2 + 4*s + 3
Hs = Ps/Qs
Siendo H(s) la forma en que se expresan los ejercicios, se realizan las operaciones para obtener las raíces del numerador y denominador. El resultado buscado con el algoritmo es:
H(s):
2*(s + 2)
---------------
(s + 1)*(s + 3)
H(s) en fracciones parciales
1 1
----- + -----
s + 3 s + 1
{Q_polos:veces}: {-1: 1, -3: 1}
{P_ceros:veces}: {-2: 1}
>>>
Instrucciones en Python
Hs se convierte en una expresión de Sympy con sym.sympify(Hs)s para el caso en que se de tan solo una constante. En el caso que Hs requiera agrupar términos se usa la instrucción sym.simplify(Hs,inverse=True).
Una vez preparada la expresión Hs, se separa como numerador Ps y denominador Qs en forma de polinomio para obtener las raíces y las veces que ocurren con la instrucción sym.roots().
La expresión de H(s) en fracciones parciales se obtienen usando sym.apart(Hs,s).
# Busca Polos y Ceros de H(s) con Sympy# Ps es numerador, Qs es denominadorimport sympy as sym
# INGRESO
s = sym.Symbol('s')
Ps = 2*s + 4
Qs = s**2 + 4*s + 3
Hs = Ps/Qs
# PROCEDIMIENTO
Hs = sym.sympify(Hs) # convierte a sympy una constante
Hs = sym.simplify(Hs,inverse=True) # agrupa terminos# polos y ceros de Hs
[P,Q] = Hs.as_numer_denom()
P = P.as_poly(s) # numerador
Q = Q.as_poly(s) # denominador
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)
# en factores
Hs = sym.factor(Hs,s)
# fracciones parciales
Hs_fp = sym.apart(Hs,s)
# SALIDAprint('H(s):')
sym.pprint(Hs)
print('\n H(s) en fracciones parciales ')
sym.pprint(Hs_fp)
print('\n {Q_polos:veces}:',Q_polos)
print('\n {P_ceros:veces}:',P_ceros)
Para visualizar el concepto de polos y ceros, se presenta la gráfica de polos y ceros en el dominio s, para un corte en el plano real. Por simplicidad se usará el procedimiento creado para el curso en telg1001.py
Todos los polos y ceros se encuentran en el lado izquierdo del plano, es decir su parte real es negativa.
# GRAFICA -----------import matplotlib.pyplot as plt
import telg1001 as fcnm
fig_Hs = fcnm.graficar_Fs(Hs,Q_polos,P_ceros,f_nombre='H')
plt.show()
Para revisar los polos de la expresión de fracciones parciales, se considera que los términos de las expresiones podrían tener polos repetidos como en el ejemplo 3. Para polos repetidos, se usa el número de veces mayor cuando aparecen polos de términos suma.
Al algoritmo del ejercicio anterior se modifica solo el bloque de ingreso con las expresiones para numerador y denominador:
Ps = 2*s**2+5
Qs = s**2 + 3*s + 2
Hs = Ps/Qs
El resultado Hs_fp en fracciones parciales muestra que existe un término constante en las sumas de la expresión. La forma descrita por Lathi al inicio de la página, la constante se la considera como ‘ganancia(s)’.
H(s):
2
2*s + 5
---------------
(s + 1)*(s + 2)
H(s) en fracciones parciales
13 7
2 - ----- + -----
s + 2 s + 1
{Q_polos:veces}: {-1: 1, -2: 1}
{P_ceros:veces}: {-sqrt(10)*I/2: 1, sqrt(10)*I/2: 1}
>>>
Por otra parte se tiene que los ceros tienen raíces con componente imaginario representado con el símbolo ‘I‘, que para representar cada raíz en la gráfica se considera el eje vertical también como parte Imag(s). Si es cierto, es diferente a la naturaleza de H(s), perdonen, se usa ésta «libertad» solo para resaltar la influencia de los polos en H(s) como se mostró en el ejercicio anterior, o se realizaría un gráfico con el plano real e imaginario como base y una tercera dimensión para H(s).
Los polos se encuentran en el lado izquierdo del plano, es decir su parte real es negativa.
Ejemplo 3. Polos reales repetidos o raices repetidas en denominador Q(s)
Referencia: Lathi ejemplo 4.3d p342
H(s)=(s+1)(s+2)38s+10
la forma de la expresión para el algoritmo es
Ps = 8*s+10
Qs = (s+1)*((s+2)**3)
Hs = Ps/Qs
En él ejercicio el término del denominador tiene un factor elevado al cubo, lo que muestra que la raiz es repetida y se debe reflejar en la respuesta del algoritmo. En las fracciones parciales se muestra que para el polo repetido existe un término suma para cada grado en el denominador.
H(s):
2*(4*s + 5)
----------------
3
(s + 1)*(s + 2)
H(s) en fracciones parciales
2 2 6 2
- ----- - -------- + -------- + -----
s + 2 2 3 s + 1
(s + 2) (s + 2)
{Q_polos:veces}: {-1: 1, -2: 3}
{P_ceros:veces}: {-5/4: 1}
>>>
Todos los polos y ceros se encuentran en el lado izquierdo del plano s.
Se adjunta la gráfica para la interpretación de H(s)
Si la busqueda de polos es sobre la expresión de fracciones parciales, se tiene que las veces se determinan como el mayor exponente del término del denominador.
Ejemplo 4. H(s) con Fracciones parciales y polos de tipo complejo
Referencia: Lathi ejemplo 4.3c p340
H(s)=s(s2+10s+34)6(s+34)
la forma de la expresión para el algoritmo es
Ps = 6*(s+34)
Qs = s*(s**2+10*s+34)
Hs = Ps/Qs
Incorporando las instrucciones para identificar si el denominador tiene raices complejas, se asignan los valores de los parametros necesarios en un diccionario de variables, obteniendo como resultado,
Ejemplo 5. Con suma de términos y exponenciales de retraso en tiempo
H(s)=s−1s+2+s−2e−2s
la expresión para el algoritmo se escribe como:
Hs = ((s+2)/(s-1)) + sym.exp(-2*s)/(s-2)
Para el segundo término que tiene un exponencial, no se puede aplicar directamente la instrucción sym.apart(Hs,s) al no ser reconocida como tipo polinomio. Para continuar se debe separar el término sym.exp() del resto de la expresión, y se podrá aplicar fracciones parciales a esa parte, para luego volver a poner el término.
Para facilitar el desarrollo del algoritmo, se recomienda empezar a analizar una expresión mas simple, como solamente el segundo término de la suma.
expresion H(s):
-2*s
s + 2 e
----- + -----
s - 1 s - 2
Polos y ceros:
término: 1 * (s + 2)/(s - 1)
Q_polos{polos:veces}: {1: 1}
P_ceros{polos:veces}: {-2: 1}
término: exp(-2*s) * 1/(s - 2)
Q_polos{polos:veces}: {2: 1}
P_ceros{polos:veces}: {}
>>>
En el algoritmo se empieza por determinar si la expresión de Hs tiene términos sym.exp(-s), la instrucción que permite listar estos términos es
lista_exp = list(Hs.atoms(sym.exp(-s)))
como la expresión H(s) en los ejercicios puede contener sumas o multiplicaciones de bloques, primero se simplifica la expresión H(s)
# convierte a Sympy si es solo constante
Hs = sym.sympify(Hs)
# simplifica operación Hs1(s)*Hs2(s)
Hs = sym.expand_power_exp(Hs)
Hs = sym.expand(Hs, power_exp=False)
term_sum = sym.Add.make_args(Hs)
Hs_0 = 0*s # reagrupa Fsfor term_k in term_sum:
term_k = sym.simplify(term_k)
term_k = sym.expand_power_exp(term_k)
Hs_0 = Hs_0 + term_k
# tabula sym.exp(-s) en lista_exp
lista_exp = list(Hs_0.atoms(sym.exp(-s)))
iflen(lista_exp)>0: # elimina constantesfor k in lista_exp:
ifnot(k.has(s)): # es constante
lista_exp.remove(k)
luego se tabula la lista_exp y se eliminan los términos que pueden ser constantes o no dependen de ‘s’, por ejemplo sym.exp(10).
Luego se separa la expresión por grupos usando sym.collect() creando un diccionario con los elementos de la lista_exp y con las expresiones agrupadas.
# separados por lista_exp
separados = sym.collect(Hs_0,lista_exp,evaluate=False)
Con los elementos separados se puede aplicar la búsqueda de polos y ceros para cada elemento de la lista con la función polosceros_simple(Hs_k)
Al final se pueden agrupar los Q_polos y P_ceros de toda la expresión siguiendo criterios como que las veces que ocurre cada polo es el maximo entre los términos separados por lista_exp.
Instrucciones en Python
# Busca Polos y Ceros de H(s) con Sympy# separa términos con sym.exp(-a*s)# Ps es numerador, Qs es denominadorimport numpy as np
import matplotlib.pyplot as plt
import sympy as sym
import telg1001 as fcnm
# INGRESO
s = sym.Symbol('s')
# Ps = (2*s + 4)# Qs = s**2 + 4*s + 3# Ps = 2*s**2+5# Qs = s**2 + 3*s + 2# Ps = 8*s+10# Qs = (s+1)*((s+2)**3)# Ps = 6*(s+34)# Qs = s*(s**2+10*s+34)# Hs = Ps/Qs# Hs = sym.exp(-2*s)/(s-2)# Hs = 1/s - sym.exp(-2*s)/s# Hs = sym.exp(-2*s)/s -sym.exp(-4*s)/s# Hs = (2*s**2+5*s+12)/((s+2)*(s**2+2*s+10))
Hs = ((s+2)/(s-1)) + sym.exp(-2*s)/(s-2)
# PROCEDIMIENTOdefbusca_polosceros(Hs):
''' Busca polos de Hs (divisiones para cero)
y ceros de Hs,cuenta las veces que aparecen,
agrupa por exp(-a*s) agrupado en factores.
'''defpolosceros_simple(Hs):
''' Busca_polo de un termino sin exp(-a*s)
'''
Hs = sym.simplify(Hs,inverse=True)
# Hs debería ser un solo termino suma
Hs_factor = sym.factor(Hs,s)
# polos y ceros de termino Hs
[P,Q] = Hs_factor.as_numer_denom()
P = sym.poly(P,s) # numerador
Q = sym.poly(Q,s) # denominador
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)
respuesta = {'Q_polos' : Q_polos,
'P_ceros' : P_ceros}
return(respuesta)
# convierte a Sympy si es solo constante
Hs = sym.sympify(Hs)
# simplifica operación Hs1(s)*Hs2(s)
Hs = sym.expand_power_exp(Hs)
Hs = sym.expand(Hs, power_exp=False)
term_sum = sym.Add.make_args(Hs)
Hs_0 = 0*s # reagrupa Fsfor term_k in term_sum:
term_k = sym.simplify(term_k)
term_k = sym.expand_power_exp(term_k)
Hs_0 = Hs_0 + term_k
# tabula sym.exp(-s) en lista_exp
lista_exp = list(Hs_0.atoms(sym.exp(-s)))
iflen(lista_exp)>0: # elimina constantesfor k in lista_exp:
ifnot(k.has(s)): # es constante
lista_exp.remove(k)
# separados por lista_exp
separados = sym.collect(Hs_0,lista_exp,evaluate=False)
# polos y ceros por terminos exp(-s) agrupados
polosceros = {} ; Hs_fp = 0
for k in separados:
Hs_k = sym.factor(separados[k],s)
PZ_k = polosceros_simple(Hs_k)
PZ_k['Hs_k'] = Hs_k
polosceros[k] = PZ_k
# integra polos
Q_polos = {} ; P_ceros = {}
for k in polosceros:
polos = polosceros[k]['Q_polos']
for unpolo in polos:
if unpolo in Q_polos:
veces = max([Q_polos[unpolo],polos[unpolo]])
Q_polos[unpolo] = veces
else:
Q_polos[unpolo] = polos[unpolo]
ceros = polosceros[k]['P_ceros']
for uncero in ceros:
if uncero in P_ceros:
veces = max([P_ceros[uncero],ceros[uncero]])
P_ceros[unpolo] = veces
else:
P_ceros[uncero] = ceros[uncero]
# revisa exp(-a*s)iflen(lista_exp)==0: #sin componentesdel polosceros['1']
polosceros['Q_polos'] = Q_polos
polosceros['P_ceros'] = P_ceros
return(polosceros)
polosceros = busca_polosceros(Hs)
# SALIDAprint('\n expresion H(s):')
sym.pprint(Hs)
print('\nPolos y ceros:')
lista_PZ = ['Q_polos','P_ceros']
for k in polosceros:
ifnot(k in lista_PZ):
print('término:', k,'*',polosceros[k]['Hs_k'])
print(' Q_polos{polos:veces}: ',polosceros[k]['Q_polos'])
print(' P_ceros{polos:veces}: ',polosceros[k]['P_ceros'])
else:
print(k,polosceros[k])
La gráfica en éste caso se realiza usando Hs y sus componentes, de ésta forma es posible visualizar el efeto de cada componente con sus polos y ceros.
# GRAFICA ----------------------defgraficar_Fs2(Fs,Q_polos={},P_ceros={},
s_a=1,s_b=0,muestras=101,f_nombre='F',
solopolos=False,polosceros={}):
# grafica Fs
lista_Hs =[]
lista_PZ = ['Q_polos','P_ceros']
lista_term = list(polosceros.keys())
lista_term.remove('Q_polos')
lista_term.remove('P_ceros')
# intervalo para s,y componentes Fs con exp(-a*s)
s_a = 0 ; s_b = 0 ; Lista_Hs = []
cond_componente = Falseiflen(lista_term)>0:
cond_componente = Trueiflen(polosceros)>1: # tiene componenes con exp(-a*s)
lista_Hs.append([Hs,{},{}])
for k in polosceros:
ifnot(k in lista_PZ):
Q_polosk = polosceros[k]['Q_polos']
P_cerosk = polosceros[k]['P_ceros']
[s_ak,s_bk] = fcnm.intervalo_s(Q_polosk,P_cerosk)
s_a = min(s_a,s_ak)
s_b = max(s_b,s_bk)
lista_Hs.append([polosceros[k]['Hs_k']*sym.sympify(k),
Q_polosk,P_cerosk])
else: #sin componenes con exp(-a*s)
lista_Hs.append([Hs,polosceros['Q_polos'],polosceros['P_ceros']])
[s_a,s_b] = fcnm.intervalo_s(polosceros['Q_polos'],
polosceros['P_ceros'])
# graficas por cada componente
fig_Fs, graf_Fs = plt.subplots()
# no presenta errores de división para cero en lambdify()
np.seterr(divide='ignore', invalid='ignore')
for k in lista_Hs:
Fs_k = k[0] ;Q_polosk = k[1] ;P_cerosk = k[2]
# convierte a sympy una constante
Fs_k = sym.sympify(Fs_k,s)
if Fs_k.has(s): # no es constante
F_s = sym.lambdify(s,Fs_k,modules=fcnm.equivalentes)
else:
F_s = lambda s: Fs_k + 0*s
s_i = np.linspace(s_a,s_b,muestras)
Fsi = F_s(s_i) # Revisar cuando s es complejo
lineaestilo = 'solid'iflen(Q_polosk)>0 and cond_componente:
lineaestilo = 'dashed'
graf_Fs.plot(s_i,Fsi,label=Fs_k,linestyle=lineaestilo)
lineacolor = plt.gca().lines[-1].get_color()
iflen(Q_polosk)>0:
polo_re = [] ; polo_im = []
for polos in Q_polosk.keys():
graf_Fs.axvline(sym.re(polos),color='red',
linestyle='dotted')
x_polo = sym.re(polos)
y_polo = sym.im(polos)
polo_re.append(x_polo)
polo_im.append(y_polo)
etiqueta = "("+str(x_polo)+','+str(y_polo)+")"
graf_Fs.annotate(etiqueta,(x_polo,y_polo))
graf_Fs.scatter(polo_re,polo_im,marker='x',
color =lineacolor,
label = Q_polosk)
iflen(P_cerosk)>0:
cero_re = [] ; cero_im = []
for cero in P_cerosk.keys():
x_cero = sym.re(cero)
y_cero = sym.im(cero)
cero_re.append(x_cero)
cero_im.append(y_cero)
etiqueta = "("+str(x_cero)+','+str(y_cero)+")"
graf_Fs.annotate(etiqueta,(x_cero,y_cero))
graf_Fs.scatter(cero_re,cero_im,marker='o',
color =lineacolor,
label = P_cerosk)
graf_Fs.axvline(0,color='gray')
graf_Fs.legend()
graf_Fs.set_xlabel('Real(s)')
graf_Fs.set_ylabel('F(s) ; Imag(s)')
graf_Fs.set_title(r'F(s) = $'+str(sym.latex(Hs))+'$')
graf_Fs.grid()
return(fig_Fs)
fig_Fs = graficar_Fs2(Hs,polosceros=polosceros)
plt.show()
Parametros cuadráticos
Para el ejercicio 4, se tienen polos con valores complejos y para usar la tabla de transformadas de Laplace se requiere obtener algunos parámetros.
o la otra foma de expresión usando la fila 12d de la tabla
Ejemplo 4 Desarrollo con Sympy-Python
Incorporando las instrucciones para identificar si el denominador tiene raices complejas, se asignan los valores de los parámetros necesarios en un diccionario de variables, obteniendo como resultado,
# Fracciones parciales de H(s) con Sympy# Ps es numerador, Qs es denominador# con Transformada Inversa de Laplaceimport numpy as np
import sympy as sym
import telg1001 as fcnm
# INGRESO
s = sym.Symbol('s')
t = sym.Symbol('t', real=True)
Ps = 6*(s+34)
Qs = s*(s**2+10*s+34)
Hs = Ps/Qs
# PROCEDIMIENTO# parametros de termino con Q cuadraticodefQ_cuad_s_parametros(Hs):
'''Si Q tiene grado=2, obtiene los parametros
para la tabla de Transformadas Inversas de Laplace
'''defQ_cuad_sin_exp(Hs):
'''para Hs sin sym.exp(-a*s)
'''# fracciones parciales NO separa Q con factores complejos
Hs_fp = sym.apart(Hs,s)
respuesta = {}
term_suma = sym.Add.make_args(Hs_fp)
for term_k in term_suma:
[P,Q] = term_k.as_numer_denom()
if sym.degree(Q) == 2 and sym.degree(P) == 1:
P_coef = P.as_coefficients_dict(s)
Q_coef = Q.as_coefficients_dict(s)
Q0 = 1 # Normalizar coefficiente s**2=1if Q_coef[s**2]!=1:
Q0 = Q_coef[s**2]
# Parametros de Q cuadratico
a = 0 ; c = 0
if s**1 in Q_coef:
a = float(Q_coef[s**1]/Q0)/2
if s**0 in Q_coef:
c = float(Q_coef[s**0]/Q0)
A = float(P_coef[s**1])
B = 0
if s**0 in P_coef:
B = float(P_coef[s**0])
rP = (A**2)*c + B**2 - 2*A*B*a
rQ = c - a**2
r = np.sqrt(rP/rQ)
b = np.sqrt(c-a**2)
thetaP = A*a-B
thetaQ = A*np.sqrt(c-a**2)
theta = np.arctan(thetaP/thetaQ)
parametro = {'A': A, 'B': B,'a': a, 'c': c,
'r': r, 'b': b,'theta':theta}
respuesta[term_k] = parametro
return (respuesta)
Hs = sym.sympify(Hs,s)
Hs_fp = fcnm.apart_s(Hs) # separa por sym.exp(-a*s)# tabula sym.exp(-s) en lista_exp
lista_exp = list(Hs_fp.atoms(sym.exp(-s)))
iflen(lista_exp)>0: # elimina constantesfor k in lista_exp:
ifnot(k.has(s)): # es constante
lista_exp.remove(k)
# separados por lista_exp
separados = sym.collect(Hs_fp,lista_exp,evaluate=False)
Qs2 = {} # agrupa parametros Qs2for k in separados:
Qs2_k = Q_cuad_sin_exp(separados[k])
iflen(Qs2_k)>0:
for term in Qs2_k:
Qs2[k*term] = Qs2_k[term]
return(Qs2)
Hs_fp = fcnm.apart_s(Hs)
Qs2 = Q_cuad_s_parametros(Hs)
# SALIDAprint('H(s) en fracciones parciales:')
sym.pprint(Hs_fp)
# parametros de termino con Q cuadraticoiflen(Qs2)>0:
print('\n parametros de Q_cuadratico: ')
for k in Qs2:
print(' ',k,' :')
print(' ',Qs2[k])
