Etiqueta: ejercicios resueltos python

Ejercicio: 2Eva2016TII_T4 resolver en dominio de frecuencia

Partiendo de,

$y(t) = x(t) \circledast h(t)$ $Y(\omega) = X(\omega) H(\omega)$ $Y\Big(\frac{\omega}{3}\Big) = X\Big(\frac{\omega}{3}\Big) H\Big(\frac{\omega}{3}\Big)$

se tiene también que:

$g(t) = x(3t) \circledast h(3t)$ $G(\omega) = \frac{1}{3}X\Big(\frac{\omega}{3}\Big)\text{ } \frac{1}{3}H\Big(\frac{\omega}{3}\Big)$ $= \frac{1}{9}X\Big(\frac{\omega}{3}\Big)H\Big(\frac{\omega}{3}\Big) = \frac{1}{9}Y\Big(\frac{\omega}{3}\Big)$ $G(\omega) = \frac{1}{3} \Bigg[\frac{1}{3}Y\Big(\frac{\omega}{3}\Big)\Bigg]$ $g(t) = \frac{1}{3} y(3t)$

A = 1/3
B = 3

 

Ejercicio: 2Eva2016TII_T3 LTI DT sistemas en serie

literal a

$S1: w[n] = \frac{1}{2}w[n-1] + x[n]$ $S2: y[n] = \alpha y[n-1] + \beta w[n]$ $SG: y[n] = -\frac{1}{8}y[n-2] + \frac{3}{4}y[n-1]+x[n]$

usando transformada z:

S1:

$W(z) = \frac{1}{2} z^{-1} W(z) + X(z)$ $W(z) - \frac{1}{2} z^{-1} W(z) = X(z)$ $W(z)\Big[1 - \frac{1}{2} z^{-1} \Big] = X(z)$ $W(z) = \frac{ X(z)}{\Big[1 - \frac{1}{2} z^{-1} \Big]}$

S2:

$Y(z) = \alpha z^{-1} Y(z) + \beta W(z)$ $Y(z) - \alpha z^{-1} Y(z) = \beta W(z)$ $Y(z) \Big[1 - \alpha z^{-1} \Big] = \beta W(z)$

sustituyendo la ecuacion de S1 para W(z)

$Y(z) \Big[1 - \alpha z^{-1} \Big] = \beta \frac{ X(z)}{\Big[1 - \frac{1}{2} z^{-1} \Big]}$ $Y(z) \frac{1}{\beta} \Big[1 - \alpha z^{-1} \Big]\Big[1 - \frac{1}{2} z^{-1} \Big] = X(z)$ $Y(z) \Big[ \frac{1}{\beta} - \frac{1/2+\alpha}{\beta} z^{-1}+ \frac{1}{2} \frac{\alpha}{\beta} z^{-2}\Big]= X(z)$

SG:

$Y(z) = -\frac{1}{8} z^{-2} Y(z) + \frac{3}{4} z^{-1}Y(z)+X(z)$ $Y(z)+\frac{1}{8} z^{-2} Y(z) - \frac{3}{4} z^{-1}Y(z) = X(z)$ $Y(z)\big[ 1 - \frac{3}{4} z^{-1} +\frac{1}{8} z^{-2} \Big] = X(z)$

comparando con la ecuación de S2

$\frac{1}{\beta} = 1$ $\beta = 1$ $\frac{1}{2} \frac{\alpha}{\beta} = \frac{1}{8}$ $\alpha = \frac{2\beta}{8} = \frac{1}{4}$

comprobar con

$- \frac{1/2+\alpha}{\beta} = \frac{3}{4}$

se confirma que  α = 1/4 y β=1

la función de transferencia es:

$Y(z)\Big[ 1 - \frac{3}{4} z^{-1} +\frac{1}{8} z^{-2} \Big] = X(z)$ $\frac{Y(z)}{X(z)} = \frac{1}{\Big[ 1 - \frac{3}{4} z^{-1} +\frac{1}{8} z^{-2} \Big]}$ $H(z) = \frac{z^2}{z^2 - \frac{3}{4} z +\frac{1}{8}}$ $\frac{H(z)}{z} = \frac{z}{z^2 - \frac{3}{4} z +\frac{1}{8}}$

usando las raices para:

$z^2 - \frac{3}{4} z +\frac{1}{8} = \Big[ z-\frac{1}{4}\Big] \Big[z-\frac{1}{2}\Big]$

y la parte derecha de la ecuación:

$\frac{z}{\Big[ z-\frac{1}{4}\Big] \Big[z-\frac{1}{2}\Big]} = \frac{C_1}{z-\frac{1}{4}}+\frac{C_2}{z-\frac{1}{2}}$

despejando para C1 y haciendo z=1/4,

$C_1 = \frac{z}{\Big[z-\frac{1}{2}\Big]} = \frac{\frac{1}{4}}{\Big[\frac{1}{4}-\frac{1}{2}\Big]} = -1$

despejando para C2 y haciendo z=1/2,

$C_2 = \frac{z}{\Big[z-\frac{1}{4}\Big]} = \frac{\frac{1}{2}}{\Big[\frac{1}{2}-\frac{1}{4}\Big]} = 2$

se H(z) se resume en,

$H(z) = - \frac{z}{z - \frac{1}{4}} +2\frac{z}{z-\frac{1}{2}}$

Para obtener h[n] usando la antitransformada,

$h[n]=z^{-1} \Big[ H(z) \Big]$ $h[n]=z^{-1} \Big[- \frac{z}{z - \frac{1}{4}} +2\frac{z}{z-\frac{1}{2}} \Big]$ $h[n]= 2 \Big[ \frac{1}{2}\Big]^n \mu [n] + \Big[ \frac{1}{4}\Big]^n \mu [n]$ $h[n]= \Bigg[ 2 \Big[ \frac{1}{2}\Big]^n - \Big[ \frac{1}{4}\Big]^n \Bigg] \mu [n]$

siendo la forma de la respuesta un impulso, es un sistema IIR.

---

Algoritmo en Python

Usando la expresión H(z) se obtiene:

 Hz:
      2     
     z      
------------
 2   3*z   1
z  - --- + -
      4    8

 Hz en fracciones parciales
     z        2*z  
- ------- + -------
  z - 1/4   z - 1/2

 Hz en factores
          2         
         z          
--------------------
(z - 0.5)*(z - 0.25)

 {Q_polos:veces}: {1/2: 1, 1/4: 1}
 {P_ceros:veces}: {0: 2}

estabilidad asintótica en z:
circ1_dentro : 2
circ1_repetidos : 0
circ1_sobre : 0
circ1_fuera : 0
unicos : 2
repetidos : 0
asintota : estable

 h[n]:
/      n        n\             
\- 0.25  + 2*0.5 /*Heaviside(n)
>>> 

añadiendo instrucciones para graficar h[n] se obtiene

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]
h[n]: [1.         0.75       0.4375     0.234375   0.12109375
       0.06152344 0.03100586 0.01556396 0.00779724 0.00390244]

Instrucciones en Python

strong>Nota: cuando se produzca el siguiente error con Numpy para evaluar una expresión con exponente negativo,

Traceback (most recent call last):
  File "D:\MATG1052Ejemplos\Transformadaz\ejercicio....py", line 93, in 
    fi  = f_n(ki)
  File "", line 2, in _lambdifygenerated
    return (9/59)*4**(-n)*Heaviside(n, 1/2)
ValueError: Integers to negative integer powers are not allowed.

proceda actualizando los valores a evaluar como tipo real (dtype float), tan solo usando en la línea de ki con lo siguiente:

ki  = np.arange(0,muestras_fn,1.0)

quedando las instrucciones de la siguiente forma, que si evalua valores para realizar gráficas.

# Transformada z- Fracciones parciales
# https://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(3,4)
a1 = sym.Rational(1,8)

Pz = z**2
Qz = z**2-a0*z+a1

#Pz = z*z**2
#Qz = (z-1)*(z**2-(a0)*z+a1)

F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())
Fz_factor = fcnm._round_float_is_int(Fz_factor)

# polos y ceros de Hz
[P,Q] = Fz.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)

estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; Fz_revisar = [] ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))
fn = fcnm._round_float_is_int(fn)

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if len(Fz_revisar)>0:
    print('revisar terminos sin transformada de tabla:')
    for un_term in Fz_revisar:
        print(un_term)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

 

---

Literal c

revisando los polos y ceros:

ceros : z = 0  y z = 0
polos: z=1/4 y z=1/2

Dado que todos los polos se encuentran dentro del círculo de radio unitario , el sistema es asintóticamente estable, por lo que es BIBO o EASA estable

literal d

$S(z) = \Big[ \frac{z}{z-1} \Big] \frac{z^2}{z^2-\frac{3}{4}z + \frac{1}{8}}$ $S(z) = \frac{z^3}{(z-1)(z-1/4)(z-1/2)}$

aplicando el mismo método anterior, se tiene:

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(3.4)
a1 = sym.Rational(1,8)

Pz = z*z**2
Qz = (z-1)*(z**2-(a0)*z+a1)
F = Pz/Qz

# para graficar
f_nombre = 'S'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

con resultado:

 Sz:
           3          
          z           
----------------------
        / 2   3*z   1\
(z - 1)*|z  - --- + -|
        \      4    8/

 Sz en fracciones parciales
     z          2*z        8*z   
----------- - ------- + ---------
3*(z - 1/4)   z - 1/2   3*(z - 1)

 Sz en factores
              3             
             z              
----------------------------
(z - 1)*(z - 0.5)*(z - 0.25)

 {Q_polos:veces}: {1: 1, 1/2: 1, 1/4: 1}
 {P_ceros:veces}: {0: 3}

 s[n]:
/    n             \             
|0.25         n   8|             
|----- - 2*0.5  + -|*Heaviside(n)
\  3              3/             

$S(z) = \frac{1}{3}\frac{z}{z-1/4} -2\frac{z}{z-1/2} + \frac{8}{3}\frac{z}{z-1}$

aplicando la transformada inversa

$s(n)= \Bigg[ \frac{1}{3}\Big[ \frac{1}{4} \Big]^n -2\Big[\frac{1}{2} \Big]^n + \frac{8}{3} \Bigg] \mu [n]$

señal discreta s[n]
n   : [0 1 2 3 4 5 6 7 8 9]
h[n]: [1.         1.75       2.1875     2.421875   
 2.54296875 2.60449219 2.63549805 2.65106201 2.65885925
 2.66276169]

 

Ejercicio: 2Eva2016TII_T2 LTI CT Circuito RC respuesta de frecuencia H(ω), impulso h(t)

literal a

$v_1 (t) = v_R (t) +v_C (t)$ $v_1 (t) = R i(t) +v_2 (t)$ $i(t) = C \frac{\delta v_2 (t)}{\delta t}$ $v_2(t) = RC \frac{\delta v_2 (t)}{\delta t} +v_2 (t)$ $V_1 (\omega) = j \omega RC V_2 (\omega) + V_2(\omega) = V_2(\omega) [1+j\omega RC]$ $H(\omega) = \frac{V_2 (\omega)}{V_1(\omega)} = \frac{1}{1+j\omega RC} = \frac{1}{1+j\frac{\omega}{\omega_c}}$ $\omega _c = \frac{1}{RC}$ $H(\omega) =\frac{V_2(\omega)}{V_1(\omega)} = \frac{1}{1+j\frac{\omega}{2}}$ $\begin{cases} |H(\omega) = \frac{1}{\sqrt{1+\big( \frac{\omega}{2}\big)^2}} \\ \theta_{H(\omega)} = -tg^{-1} \big( \frac{\omega}{2}\big) \end{cases}$ $\omega _c = \frac{1}{RC}$

literal c

La respuesta impulso del filtro LPF, se obtiene mediante:

$h(t) = \mathcal{F}^{-1} \Big[ H(\omega) \Big] = \mathcal{F}^{-1} \Big[ \frac{1}{1+j(\omega /2)} \Big]$ $=2\mathcal{F}^{-1} \Big[ \frac{1}{j\omega +2} \Big]$ $h(t)=2e^{2t}\mu (t)$

literal d

Método 1: usando la respuesta de frecuencia

$\begin{cases} |H(\omega) = \frac{1}{\sqrt{1+\big( \frac{\omega}{2}\big)^2}} \\ \theta_{H(\omega)} = -tg^{-1} \big( \frac{\omega}{2}\big) \end{cases}$ $\begin{cases} |H(50) = \frac{1}{\sqrt{1+\big( \frac{50}{2}\big)^2}} = \frac{1}{\sqrt(626)} \\ \theta_{H(\omega)} = -tg^{-1} \big( \frac{50}{2}\big) = -87.7\end{cases}$ $v_2(t) = |H(50)| \sin \big(50t+\theta_{H(50)} \big)$ $v_2(t) = \frac{1}{\sqrt(626)} \sin \big(50t-87.7) \big)$

método 2:

$V_2(\omega ) = V_1(\omega) H(\omega)$ $V_1 (\omega) = \mathcal{F}[v_1(t) ] = \mathcal{F} [ \sin (50t) ]$ $= j \pi \delta (\omega +50) - j \pi \delta (\omega-50)$ $V_2 (\omega) = V_1(\omega) H(\omega)$ $V_2 (\omega) = \Big[ j \pi \delta (\omega +50) - j \pi \delta (\omega-50) \Big] \Big[ \frac{1}{1+j(\omega/2)} \Big]$ $= j \pi \Big[ \delta (\omega + 50)\frac{1}{1+j(\omega /2)} - \delta (\omega - 50)\frac{1}{1+j(\omega /2)}\Big]$ $= j \pi \Big[ \delta (\omega + 50)\frac{1}{1+j(50 /2)} - \delta (\omega - 50)\frac{1}{1+j(50 /2)}\Big]$ $= j \pi \Big[ \delta (\omega + 50)\frac{1+j25}{626} - \delta (\omega - 50)\frac{1-j25}{626}\Big]$ $= j \frac{\pi}{626} \Big[ \delta (\omega + 50) - \delta (\omega - 50) + j25 \delta (\omega + 50) + j25 \delta (\omega - 50)\Big]$ $= \frac{1}{626} \Big[ j\pi \delta (\omega + 50) - j \pi \delta(\omega - 50)\Big] - \frac{25}{626} \Big[ \pi \delta (\omega + 50) +\pi \delta(\omega - 50)\Big]$ $v_2(t) = \mathcal{F}^{-1} [V_2 (\omega)]$ $= \frac{1}{626} \mathcal{F}^{-1}\Bigg[ \Big[ j\pi \delta (\omega + 50) - j \pi \delta(\omega - 50)\Big] - 25 \Big[ \pi \delta (\omega + 50) +\pi \delta(\omega - 50)\Big] \Bigg]$ $v_2(t) = \frac{1}{626} \Big[ \sin (50t) - 25 \cos (50t) \Big]$

usando fasores:

$v_2(t) = \frac{1}{626} \cos (50t-177.709) = \frac{1}{626} \sin (50t-87.70)$

En la salida, existe un factor de atenuación de 0.04 y un retardo de 87.70°.
Como la señal de entrada se reproduce de manera exacta en su salida a pesar tener amplitud diferente y un retardo en el tiempo.

Ejercicio: 2Eva2012TI_T3 LTI DT causal, coeficientes de respuesta impulso h[n]

La ecuación de diferencias descrita se usa para crear el diagrama de bloques.

$y[n] - \frac{5}{4} y[n-1] + \frac{1}{36} y[n-2] + \frac{1}{18} y[n-3] = x[n] - \frac{1}{2} x[n-1]$

literal a. Respuesta impulso h[n]

La forma para usar transformadas z requiere que se desplazar tres unidades.

$y[n+3] - \frac{5}{4} y[n+2] + \frac{1}{36} y[n+1] + \frac{1}{18} y[n] = x[n+3] - \frac{1}{2} x[n+2]$ $z^3 Y[z] - \frac{5}{4} z^2 Y[z] + \frac{1}{36} z Y[z] + \frac{1}{18} Y[z] = z^3X[z] - \frac{1}{2} z^2 X[z]$ $(z^3 - \frac{5}{4} z^2 + \frac{1}{36} z + \frac{1}{18}) Y[z] = (z^3 - \frac{1}{2} z^2) X[z]$ $H[z] = \frac{Y[z]}{X[z]} = \frac{z^3 - \frac{1}{2} z^2}{z^3 - \frac{5}{4} z^2 + \frac{1}{36} z + \frac{1}{18}}$

las raíces del denominador obtenidas con Sympy-Python:

>>> sym.roots(Pz)
{1/2: 1, 0: 2}
>>> sym.factor(Pz)
z**2*(2*z - 1)/2
>>> sym.factor(Qz)
(4*z - 1)*(9*z**2 - 9*z - 2)/36
>>> sym.roots(Qz)
{1/4: 1, 1/2 - sqrt(17)/6: 1, 1/2 + sqrt(17)/6: 1}
>>> 

Observación: existe una raíz con distancia al origen superior al radio=1, por lo que habría un termino creciente no acotado.

que se observa también como:

El modelo se puede plantear como:

$H[z] = \frac{z^2\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)}$

para usar fracciones parciales modificadas, se multiplica cada lado por 1/z:

$\frac{H[z]}{z} = \frac{z\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)}$

Se puede plantear un modelo de respuesta por cada raíz como:

$\frac{H[z]}{z} = \frac{k1}{z-\frac{1}{4}}+\frac{k_2 z +k_3}{z^2-z-\frac{2}{9}}$

usando el método de Heaviside,

$k_1 = \frac{z\Big(z-\frac{1}{2}\Big)}{\cancel{\Big(z-\frac{1}{4}\Big)}\Big(z^2-z-\frac{2}{9}\Big)}\Big|_{z=\frac{1}{4}}$ $k_1 = \frac{\frac{1}{4}\Big(\frac{1}{4}-\frac{1}{2}\Big)}{\Big(\frac{1}{4}^2-\frac{1}{4}-\frac{2}{9}\Big)}$ $k_1 = -\frac{1}{16}\frac{1}{\frac{(4)(9)-(16)(9)-(2)(16)(4)}{(16)(4)(9)}} =\frac{36}{236}$

con lo que ahora en la expresión se convierte en:

$\frac{H[z]}{z} = \frac{z\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)} =\frac{\frac{36}{236}}{z-\frac{1}{4}}+\frac{k_2 z +K_3}{z^2-z-\frac{2}{9}} =$

Usando el método de los factores cuadráticos, se multiplica ambos lados por z y z→∞

$\frac{H[z]}{z} z = \frac{z^2\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)} =\frac{\frac{36}{236}z}{z-\frac{1}{4}}+z\frac{k_2 z +K_3}{z^2-z-\frac{2}{9}} =$

se divide numerador y denominador del lado izquierdo para 1/z³ y el lado derecho el primer termino 1/z y el segundo termino 1/z²

$\frac{\Big(1-\frac{1}{2}\frac{1}{z}\Big)}{\Big(1-\frac{1}{4}\frac{1}{z}\Big)\Big(1-\frac{1}{z}-\frac{2}{9}\frac{1}{z^2}\Big)} =\frac{\frac{36}{236}}{1-\frac{1}{4}\frac{1}{z}}+\frac{k_2 +K_3\frac{1}{z}}{1-\frac{1}{z}-\frac{2}{9}\frac{1}{z}}$

y cuando z→∞

$\frac{\Big(1-\frac{1}{2}(0)\Big)}{\Big(1-\frac{1}{4}(0)\Big)\Big(1-(0)-\frac{2}{9}(0)\Big)} =\frac{\frac{36}{236}}{1-\frac{1}{4}(0)}+\frac{k_2 +K_3(0)}{1-(0)-\frac{2}{9}(0)}$ $1 =\frac{36}{236}+k_2$ $k_2 = 1-\frac{36}{236}=\frac{200}{236} = \frac{50}{59}$

con lo que K₂=50/59 , para encontrar K₃ se usa un valor conveniente de z=0

$H[0] = \frac{(0)\Big((0)-\frac{1}{2}\Big)}{\Big((0)-\frac{1}{4}\Big)\Big((0)^2-(0)-\frac{2}{9}\Big)} =\frac{\frac{36}{236}}{(0)-\frac{1}{4}}+\frac{\frac{50}{59}(0) +K_3}{(0)^2-(0)-\frac{2}{9}} =$ $0 =\frac{\frac{36}{236}}{-\frac{1}{4}}+\frac{K_3}{-\frac{2}{9}} =$ $K_3 = \frac{2}{9}\frac{\frac{36}{236}}{-\frac{1}{4}} = -4\frac{2}{9}\frac{36}{236} = -\frac{8}{59}$ $\frac{H[z]}{z} = \frac{36}{236}\frac{1}{z-\frac{1}{4}}+\frac{50/59 z -8/59}{z^2-z-\frac{2}{9}} =$ $\frac{H[z]}{z} = \frac{36}{236}\frac{1}{z-\frac{1}{4}}+\frac{2}{59}\frac{25 z -4}{z^2-z-\frac{2}{9}}$

y restaurando en fracciones parciales al multiplicar por z cada lado

$H[z] = \frac{9}{59}\frac{z}{z-\frac{1}{4}}+\frac{2}{59}\frac{z(25 z -4)}{z^2-z-\frac{2}{9}}$

los valores no se ajustan al modelo planteado en el enunciado, y existe un polo con radio>1, por lo que el sistema es creciente,y otro polo en posición r<0.

$h[n] = a \alpha^n \mu [n] + b \beta^n \mu[n] + c \gamma^n \mu [n]$

obtenga entonces los valores pertinentes.

a= 9/59 = 0.1525	b= … creciente sin cota, con polo mayor a 1	c=… con polo negativo, lo que es un término con signo alternado.
α= 1/4 = 0.25	β = 1/2 + np.sqrt(17)/6 = 1.1871842709362768	γ = 1/2 – np.sqrt(17)/6 = -0.18718427093627676

Resultados iniciales con el algoritmo

Se encuentra que hay un término que no se puede usar para encontrar la transformada z.

 Hz:
            2      
       3   z       
      z  - --      
           2       
-------------------
        2          
 3   5*z    z    1 
z  - ---- + -- + --
      4     36   18

 Hz en fracciones parciales
50*z*(z - 4/25)       9*z     
--------------- + ------------
   / 2       2\   59*(z - 1/4)
59*|z  - z - -|               
   \         9/               

 Hz en factores
               2                       
              z *(z - 0.5)             
---------------------------------------
           / 2                        \
(z - 0.25)*\z  - z - 0.222222222222222/

 {Q_polos:veces}: {1/4: 1, 1/2 - sqrt(17)/6: 1, 1/2 + sqrt(17)/6: 1}
 {P_ceros:veces}: {1/2: 1, 0: 2}

parametros cuadraticos: 
 termino: 50*z*(z - 4/25)/(59*(z**2 - z - 2/9))
r : 120.72788283210394
gamma : None
beta : None
theta : None

estabilidad asintótica en z:
circ1_dentro : 2
circ1_repetidos : 0
circ1_sobre : 0
circ1_fuera : 1
unicos : 3
repetidos : 0
asintota : inestable

 h[n]:
   -n             
9*4  *Heaviside(n)
------------------
        59        
revisar terminos sin transformada de tabla:
50*z*(z - 4/25)/(59*(z**2 - z - 2/9))

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]
h[n]: [1.52542373e-01 3.81355932e-02 9.53389831e-03 2.38347458e-03
 5.95868644e-04 1.48967161e-04 3.72417903e-05 9.31044756e-06
 2.32761189e-06 5.81902973e-07]

Instrucciones en Python

Desarrollo con las expresiones iniciales, sin corrección sobre los parámetros que se indican en el enunciado.

Nota: cuando se produzca el siguiente error con Numpy para evaluar una expresión con exponente negativo,

Traceback (most recent call last):
  File "D:\MATG1052Ejemplos\Transformadaz\ejercicio03.py", line 93, in 
    fi  = f_n(ki)
  File "", line 2, in _lambdifygenerated
    return (9/59)*4**(-n)*Heaviside(n, 1/2)
ValueError: Integers to negative integer powers are not allowed.

proceda actualizando los valores a evaluar como tipo real (dtype float), tan solo usando en la línea de ki con lo siguiente:

ki  = np.arange(0,muestras_fn,1.0)

quedando las instrucciones de la siguiente forma, que si evalúa valores para realizar gráficas.

# Transformada z- Fracciones parciales
# https://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(1,2)
a1 = sym.Rational(5,4)
a2 = sym.Rational(1,36)
a3 = sym.Rational(1,18)

Pz = z**3-a0*z**2
Qz = z**3-a1*z**2+a2*z+a3

F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())
Fz_factor = fcnm._round_float_is_int(Fz_factor)

# polos y ceros de Hz
[P,Q] = F.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)

estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; f_noeval = 0*n ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if not f_noeval==sym.S.Zero:
    print('revisar terminos sin transformada de tabla:')
    print(f_noeval)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

 

Ejercicio: 2Eva2011TI_T2 LTI DT Determinar H[z] desde bloques

literal a. expresión de la función de transferencia

El diagrama de bloques del enunciado se reordena de la siguiente forma:

El nuevo diagrama muestra que el sistema tiene dos sub-componentes en paralelo.

$H[z] = -\frac{(11/2)z +7}{z^2-1z-2} +\frac{z}{2} -\frac{9}{2}$

como existen varios componentes, se pueden tratar de forma separada.

$H[z] = H_1[z] +\frac{z}{2} -\frac{9}{2}$ $H_1[z] = -\frac{(11/2)z +7}{z^2-z-2}$

los polos y ceros de H₁

 {Q_polos:veces}: {2: 1, -1: 1}
 {P_ceros:veces}: {-14/11: 1}

$H_1[z] = -\frac{(11/2)z +7}{(z-2)(z+1)}$

usando fracciones parciales modificadas

$\frac{H_1[z]}{z} = -\frac{(11/2)z +7}{z(z-2)(z+1)} = \frac{k_1}{z}+\frac{k_2}{z-2}+\frac{k_3}{z+1}$

usando el algoritmo se tiene:

$\frac{H_1[z]}{z} = \frac{\frac{7}{2}}{z}-\frac{3}{z-2}-\frac{\frac{1}{2}}{z+1}$

restaurando a fracciones parciales, al multiplicar cada lado por z, se obtiene la función de transferencia H[z] completa como componentes en paralelo junto con H₁

$H_1[z] = \frac{7}{2}-\frac{3z}{z-2}-\frac{1}{2}\frac{z}{z+1}$

usando la tabla de transformadas z

$h_1[n] = \frac{7}{2}\delta[n]-3(2)^n \mu[n] -\frac{1}{2}(-1)^n \mu[n]$

Usando el algoritmo de la sección X[z] Fracciones parciales modificadas con Python para H1[z]

con entrada:

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(11,2)

Pz = -(a0*z+7)
Qz = z**2-z-2

se obtiene:

 Hz:
  11*z    
- ---- - 7
   2      
----------
 2        
z  - z - 2

 Hz en fracciones parciales
      z        3*z    7
- --------- - ----- + -
  2*(z + 1)   z - 2   2

 Hz en factores
-3.5*(0.785714285714286*z + 1) 
-------------------------------
      (0.5*z - 1)*(z + 1)      

 {Q_polos:veces}: {2: 1, -1: 1}
 {P_ceros:veces}: {-14/11: 1}

estabilidad asintótica en z:
circ1_dentro : 0
circ1_repetidos : 0
circ1_sobre : 1
circ1_fuera : 1
unicos : 2
repetidos : 0
asintota : inestable

 h[n]:
/      n       \                               
|  (-1)       n|                7*DiracDelta(n)
|- ----- - 3*2 |*Heaviside(n) + ---------------
\    2         /                       2       

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]
h[n]: [    0.     -5.5   -12.5   -23.5   -48.5
         -95.5  -192.5  -383.5  -768.5 -1535.5]

resultado que se completan con los términos de los otros componentes, al incorporar la expresión con los elementos en paralelo como se desarrolla en el siguiente literal.

Instrucciones en Python

# Transformada z- Fracciones parciales
# https://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(11,2).limit_denominator(100)

Pz = -(a0*z+7)
Qz = z**2-z-2
F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())

# polos y ceros de Hz
[P,Q] = F.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)

estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; Fz_revisar = [] ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))
#fn = fcnm._round_float_is_int(fn)

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if len(Fz_revisar)>0:
    print('revisar terminos sin transformada de tabla:')
    for un_term in Fz_revisar:
        print(un_term)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

literal b. Ecuación de diferencias

Se realiza la conversión por la suma de cada componente (en paralelo):

Para y₁[n]:

$Y_1[z] = H_1[z]X[z] = \Big[-\frac{(11/2)z +7}{z^2-z-2}\Big] X[z]$ $Y_1[z] [z^2-z-2]= -(11/2)zX[z] -7X[z]$ $y_![n+2] - y_1[n+1] -2y_1[n] = -\frac{11}{2} x[n+1] -7x[n]$

Para y₂[n]:

$Y_2[z] = H_2[z]X[z] = \frac{z}{2}X[z]$ $y_2[n] = \frac{1}{2}x[n+1]$

Para y₃[n]:

$Y_3[z] = H_3[z]X[z] = -\frac{9}{2} X[z]$ $y_3[n] = -\frac{9}{2} x[n]$

se suman las expresiones obtenidas de Y₁[z] + Y₂[z]+Y₃[z]

$y[n+2] - y[n+1] -\cancel{2y[n]} +\cancel{y[n]} + \cancel{y[n]}=$ $-\frac{11}{2} x[n+1] -7x[n] +\frac{1}{2}x[n+1]-\frac{9}{2} x[n]$

la ecuación de diferencias simplificada es;

$y[n+2] - y[n+1] = -5 x[n+1] -\frac{23}{2} x[n]$

El sistema global se puede reescribir nuevamente en z como

$z^2Y[z] - zY[z] = -5 zX[z] -\frac{23}{2} X[z]$ $z(z-1)Y[z] = (-5 z -\frac{23}{2}) X[z]$

se tiene que el sistema tienen polos en 0 y 1, que se encuentran dentro del radio 1 del plano z.

literal d. Respuesta al impulso H[z]

H[z] se obtiene a partir de la última ecuación

$H[z] = \frac{Y[z]}{H[z]} = -\frac{5 z +\frac{23}{2}}{z(z-1)}$

aplicando fracciones parciales modificadas:

$\frac{H[z]}{z} = -\frac{5 z +\frac{23}{2}}{z^2(z-1)}$ $\frac{H[z]}{z} = \frac{33}{2} \frac{1}{z} +\frac{23}{2}\frac{1}{z^2} -\frac{33}{2} \frac{1}{z-1}$

restaurando fracciones parciales al multiplicar por z

$H[z] = \frac{33}{2} +\frac{23}{2} \frac{1}{z} -\frac{33}{2} \frac{z}{z-1}$

usando la tabla de transformada z se tiene:

$h[n] = \frac{33}{2}\delta [n] +\frac{23}{2} \delta [n-1] - \frac{33}{2} \mu[n]$

usando el algoritmo con entrada:

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(11,2)
a1 = sym.Rational(9,2)
a2 = sym.Rational(23,2)
#F = -(a0*z+7)/(z**2-z-2) + z/2 -a1
F = -(5*z+a2)/(z*(z-1))

se tiene como resultado:

 Hz:
-5*z - 23/2
-----------
 z*(z - 1) 

 Hz en fracciones parciales
  16.5*z          11.5
- ------ + 16.5 + ----
  z - 1            z  

 Hz en factores
-11.5*(0.434782608695652*z + 1) 
--------------------------------
           z*(z - 1)            

 {Q_polos:veces}: {1: 1, 0: 1}
 {P_ceros:veces}: {-2.30000000000000: 1}

estabilidad asintótica en z:
circ1_dentro : 1   circ1_repetidos : 0 
circ1_sobre : 1  circ1_fuera : 0
unicos : 2  repetidos : 0
asintota : marginalmente estable

 h[n]:
16.5*DiracDelta(n) + 11.5*DiracDelta(n - 1) - 16.5*Heaviside(n)

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6.]
h[n]: [  0.   -5.  -16.5 -16.5 -16.5 -16.5 -16.5]

Ejercicio: 2Eva2010TI_T3 LTI CT simplificar H(s) por Bloques

Se ubican algunos puntos de referencia sobre el diagrama de bloques para plantear las ecuaciones, dejando para el último el bloque en serie del exponencial e^(-3s) o retraso en tiempo

$W_1(s) = X(s)$ $W_2(s) = W_1(s)-\frac{6}{s+5} W_4(s)$ $W_2(s) = X(s)- -\frac{6}{s+5} W_4(s)$ $W_3(s) = \frac{1}{s} W_2(s) = \frac{1}{s} \Big[ X(s) -\frac{6}{s+5} W_4(s) \Big]$ $W_3(s) = \frac{1}{s} X(s) -\frac{6}{s(s+5)} W_4(s)$ $W_4(s) = W_3(s) +\frac{1}{s+4}X(s)$

Se plantea la simplificación de H(s) =W4(s)/X(s), dejando el término exponencial o retraso de tiempo para el final. Por lo que se despeja W4(s)

$W_4(s) = \frac{1}{s} X(s) -\frac{6}{s(s+5)} W_4(s) +\frac{1}{s+4}X(s)$ $W_4(s) +\frac{6}{s(s+5)} W_4(s) = \frac{1}{s} X(s) +\frac{1}{s+4}X(s)$ $\Big[1 +\frac{6}{s(s+5)}\Big] W_4(s) = \Big[ \frac{1}{s} +\frac{1}{s+4} \Big] X(s)$ $\frac{s(s+5)+6}{s(s+5)} W_4(s) = \frac{(s+4)+s}{s(s+4)} X(s)$ $\frac{s^2+5s+6}{s+5} W_4(s) = \frac{2(s+2)}{s+4} X(s)$ $\frac{(s+2)(s+3)}{s+5} W_4(s) = \frac{2(s+2)}{s+4} X(s)$ $\frac{s+3}{s+5} W_4(s) = \frac{2}{s+4} X(s)$ $\frac{W_4(s)}{X(s)} = \frac{2(s+5)}{(s+4)(s+3)}$

La respuesta al impulso H₁(s) tiene polos en s=-3 y s=-4, que se encuentran en el lado izquierdo del plano imaginario s. Por lo que sus componentes en el dominio del tiempo son exponenciales decrecientes, el sistema es asintóticamente estable.

El grado del polinomio P del numerador es menor al grado del polinomio Q del numerador.

Separando en fracciones parciales para H₁(s) :

$H_1(s) =\frac{2(s+5)}{(s+4)(s+3)}$ $H_1(s) = \frac{k_1}{s+4} +\frac{k_2}{s+3} =$ $k_1 =\frac{2(s+5)}{\cancel{(s+4)}(s+3)}\Bigg|_{s=-4} = \frac{2(-4+5)}{(-4+3)} = \frac{2(1)}{-1} = -2$ $k_2 =\frac{2(s+5)}{(s+4)\cancel{(s+3)}} \Bigg|_{s=-3} = \frac{2(-3+5)}{(-3+4)} = \frac{2(2)}{1} = 4$ $H_1(s) = -\frac{2}{s+4} +\frac{4}{s+3} =$

Finalmente, añadiendo el término exponencial, que es el retraso en tiempo del sistema:

$H(s) = \Big[-\frac{2}{s+4} +\frac{4}{s+3} \Big] e^{-3s}$

se observa el comportamiento de H(s) junto a los polos en la parte real e imaginaria:

el resultado con el algoritmo para el literal a es:

 H(s) = P(s)/Q(s):
/    2       4  \  -3*s
|- ----- + -----|*e    
\  s + 4   s + 3/      
 H(s) en factores:
           -3*s
2*(s + 5)*e    
---------------
(s + 3)*(s + 4)

 h(t) :
/   9  -3*t      12  -4*t\                 
\4*e *e     - 2*e  *e    /*Heaviside(t - 3)

polosceros:
exp(-3*s) : {'Q_polos': {-3: 1, -4: 1}, 'P_ceros': {-5: 1}, 'Hs_k': 2*(s + 5)/((s + 3)*(s + 4))}
Q_polos : {-3: 1, -4: 1}
P_ceros : {-5: 1}

Estabilidad de H(s):
 n_polos_real : 2
 n_polos_imag : 0
 enRHP : 0
 unicos : 0
 repetidos : 0
 asintota : estable

---

literal b. mostrar h(t)

$H(s) = \Big[-\frac{2}{s+4} +\frac{4}{s+3} \Big] e^{-3s}$

Usando la tabla de transformadas de Laplace, para la línea 5 y la propiedad de desplazamiento en t, se obtiene:

$h(t) = \Big[-2 e^{-4(t-3)} +4e^{-3(t-3)} \Big] \mu(t-3)$

literal c. Y(s) con entrada exponencial decreciente

La señal de salida Y(s)=H(s)X(s) ante una entrada X(s)=1/(s+5) que es la transformada de Laplace de x(t)=e^-5t μ(t)

$Y(s) = H(s)x(s) = \Big[\frac{2(s+5)}{(s+4)(s+3)}\Big]\Big[\frac{1}{s+5}\Big] e^{-3s}$ $Y(s) = \frac{2}{(s+4)(s+3)}e^{-3s}$

Tarea: Realizar el desarrollo analítico y revisar los resultados del algoritmo

al realizar las fracciones parciales, se obtiene:

$Y(s) = \Big[-\frac{2}{s+4} +\frac{2}{s+3} \Big] e^{-3s}$

Con transformada inversa de Laplace, usando nuevamete la tabla de transformadas de Laplace, para la línea 5 y la propiedad de desplazamiento en t, se obtiene:

$y(t) = \Big(-2 e^{-4(t-3)} +2e^{-3(t-3)}\Big) \mu (t-3)$

El resultado del algoritmo para el literal c del ejercicio es:

 X(s): 
  1  
-----
s + 5

Respuesta entrada cero ZIR H(s) y condiciones iniciales
term_cero : 0
ZIR :
0
yt_ZIR :
0

 ZSR respuesta estado cero:
ZSR :
/    2       2  \  -3*s
|- ----- + -----|*e    
\  s + 4   s + 3/      
yt_ZSR :
/   9  -3*t      12  -4*t\                 
\2*e *e     - 2*e  *e    /*Heaviside(t - 3)

 Y(s)_total = ZIR + ZSR:
/    2       2  \  -3*s
|- ----- + -----|*e    
\  s + 4   s + 3/      

 y(t)_total = ZIR + ZSR:
/   9  -3*t      12  -4*t\                 
\2*e *e     - 2*e  *e    /*Heaviside(t - 3)
>>>

Instrucciones en 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.

# 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)

# Literal a. H(s) y estabilidad
Hs = 2*(s+5)/((s+4)*(s+3)) * sym.exp(-3*s)
#Hs = 1+0*s cuando es constante

# literal c. X(s) Señal de entrada
xt = sym.exp(-5*t)*u

# condiciones iniciales, [y'(0),y(0)] orden descendente
t0 = 0
cond_inicio = [0, 0] # estado cero no se usan

# Grafica, intervalo tiempo [t_a,t_b]
t_a = 0 ; t_b = 5
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)

# SALIDA
print(' H(s) = P(s)/Q(s):')
sym.pprint(Hs)
print(' H(s) en factores:')
sym.pprint(Hs_fc)
if len(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')

if not(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) --------
muestras_H = 201
figura_s  = fcnm.graficar_Fs(Hs_fc,Q_polos,P_ceros,f_nombre='H',solopolos=True)
figura_Hs = fcnm.graficar_Fs(Hs_fc,Q_polos,P_ceros,muestras=muestras_H,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()

Ejercicio: 2Eva2010TI_T1 LTI DT Bloques de H[z] para ecuación de diferencias

litera a. La función de transferencia H(z)

El diagrama mostrado se reordena y añade anotaciones para facilitar la escritura de la ecuación de diferencias:

Dado que hay dos bloques (1/z) o de retraso hacia abajo, se identifica que los polinomios son de segundo orden (N=2).

Los coeficientes de retraso se leen de arriba hacia abajo como a₀,-a₁,-a₂ ordenados en la siguiente tabla:

coeficientes

retraso	adelanto
a0 = 1	b0 = 0
a1 = -(1.6) = -1.6	b1 = 4
a2 = -(-0.63) = 0.63	b2 = -4

Se aplica lo mismo para los coeficientes de adelanto y se completa el modelo para los polinomios del numerador y denominador:

$H[z] = \frac{b_0 z^N +b_1 z^{N-1} + \text{ ... } + b_{N-1} z + b_N}{z^N + a_1 z^{N-1} +\text{ ... } + a_{N-1}z +a_N}$ $H[z] = \frac{4z-4}{z^2 -1.6 z + 0.63}$

Con el polinomio del denominador Qz  se encuentran las raíces

polos: {0.70: 1, 0.90: 1}

se aplican fracciones parciales modificadas al dividir cada lado para z,

$\frac{H[z]}{z} = \frac{4}{z}\frac{z-1}{(z-0.7)(z-0.9)} = \frac{k_1}{z-0.7} +\frac{k_2}{z-0.9}+\frac{k_3}{z}$

aplicando el método de Heaviside:

$k_1 = \frac{4}{z}\frac{z-1}{\cancel{(z-0.7)}(z-0.9)}\Big|_{z=0.7} = \frac{4}{0.7}\frac{(0.7)-1}{((0.7)-0.9)} = 8.5714$ $k_2 = \frac{4}{z}\frac{z-1}{(z-0.7)\cancel{(z-0.9)}}\Big|_{z=0.9} = \frac{4}{0.9}\frac{(0.9)-1}{((0.9)-0.7)} = -2.2222$ $k_3 = \frac{4}{\cancel{z}}\frac{z-1}{(z-0.7)(z-0.9)}\Big|_{z=0} = 4\frac{(0)-1}{((0)-0.7)((0)-0.9)} = -6.3492$ $\frac{H[z]}{z} = \frac{8.5714}{z-0.7} -\frac{2.2222}{z-0.9}-\frac{6.3492}{z}$

restaurando a fracciones parciales al multiplicar por z

$H[z] = \frac{8.5714z}{z-0.7} -\frac{2.2222z}{z-0.9}-6.3492$

usando la tabla de transformadas z

$h[n] = 8.5714 (0.7)^n \mu[n] -2.2222(0.9)^n \mu[n] -6.3492 \delta[n]$

Algoritmo en Python

Usando el algoritmo de la sección X[z] Fracciones parciales modificadas con Python

con entrada:

Pz = 4*z-4
Qz = z**2-1.6*z+0.63

y resultado:

 Hz:
     4*z - 4     
-----------------
 2               
z  - 1.6*z + 0.63

 Hz en fracciones parciales
8.57142857142857*z   2.22222222222222*z                   
------------------ - ------------------ - 6.34920634920635
     z - 0.7              z - 0.9                         

 Hz en factores
     4*(z - 1)     
-------------------
(z - 0.9)*(z - 0.7)

 {Q_polos:veces}: {0.700000000000000: 1, 0.900000000000000: 1}
 {P_ceros:veces}: {1: 1}

estabilidad asintótica en z:
circ1_dentro : 2
circ1_repetidos : 0
circ1_sobre : 0
circ1_fuera : 0
unicos : 2
repetidos : 0
asintota : estable

 h[n]:
/                    n                       n\                               
\8.57142857142857*0.7  - 2.22222222222222*0.9 /*Heaviside(n) 

  - 6.34920634920635*DiracDelta(n)

señal discreta h[n]
n   : [0 1 2 3 4 5 6 7 8 9]
h[n]: [ 0.   4.   2.4   1.32   0.6   0.1284
       -0.17256    -0.356988   -0.462468   -0.51504636]

Observaciones:

Las raíces características o frecuencias naturales del sistema se encuentran dentro del círculo de radio unitario. El sistema es asintóticamente estable, que implica que es BIBO estable.

literal c. Ecuación de diferencias

$H[z] = \frac{4z-4}{z^2 -1.6 z + 0.63}$ $Y[z](z^2 -1.6 z+ 0.63) = (4z-4)X[z]$ $y[n+2] - 1.6 y[n+1] + 0.63y[n] = 4x[n+1]-4x[n]$

literal d. entrada x[n]

x(t) = cos(1500t) = cos (2πf) con intervalo de muestreo T_s = 0.0015 = 1/f_s.

De la señal de entrada se tiene que 2πf=1500, por lo que f = 1500/(2π)

Para un buen muestreo de la señal de entrada, se debe cumplir la frecuencia de Nyquist, es decir f_s>2f

$fs = \frac{1}{0.0015} = 666.66 \text{ ; } f= \frac{1500}{2 \pi}=238.73$ $fs > 2*f = 2(238.73) = 477.46$

que cumple con la frecuencia mínima de muestreo y no se presentan inconvenientes para el funcionamiento de la señal por sub-muestreo

Instrucciones en Python

# Transformada z- Fracciones parciales
# https://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
Pz = 4*z-4
Qz = z**2-1.6*z+0.63

F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())
Fz_factor = fcnm._round_float_is_int(Fz_factor)

# polos y ceros de Hz
[P,Q] = F.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)
estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; Fz_revisar = [] ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))
fn = fcnm._round_float_is_int(fn)

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if len(Fz_revisar)>0:
    print('revisar terminos sin transformada de tabla:')
    for un_term in Fz_revisar:
        print(un_term)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()