```
from sympy import init_printing
init_printing()

from numpy import roots
```


```
from tabulate import tabulate
from IPython.display import Latex
from numpy import ceil, zeros

def routh_hurwitz(polinomio):
    """
    Este función toma como entrada un arreglo con los coeficientes
    del polinomio caracteristico del sistema para analizar su
    estabilidad, por medio de Routh - Hurwitz.
    """
    grado = len(polinomio)
    linea_cero = False
    long_linea = ceil((grado + 1)/2.0)
    
    def fila(tabla, i):
        """
        Se calcula la fila 'i' de una tabla de Routh - Hurwitz.
        """
        temp = []
        long_linea = ceil(grado + 1)/2.0
        for j in range(len(tabla[0]) - 1):
            # Se calcula por la formula (a1 b2 - a2 b1) / a1
            temp.append((tabla[i - 1][0]*tabla[i - 2][j + 1] -
                         tabla[i - 1][j + 1]*tabla[i - 2][0]) /
                         tabla[i - 1][0])
        # Se agrega un 0 en cada fila, por los que no se pueden calcular
        temp.append(0)
        # Se revisa el caso en que una linea se hace 0 en su totalidad
        # lo cual pasa porque la fila es linealmente dependiente con la
        # anterior.
        if (temp == zeros(long_linea)).all():
            linea_cero = True
            linea = grado - i - 1
            exponentes = [linea - j for j in range(0, linea + 1, 2)]
            while long_linea > len(exponentes):
                exponentes.append(0)
            temp =[]
            for j in range(len(tabla[0])):
                temp.append(tabla[i - 1][j]*exponentes[j])
        return temp
    
    # Se crean 2 arreglos y se introducen los coeficientes de las potencias
    # pares en uno y las impares en otro
    temp1 = []
    temp2 = []
    for i, coef in enumerate(polinomio):
        if i%2 == 0:
            temp1.append(coef)
        else:
            temp2.append(coef)
            
    # Si tiene un numero impar de potencias la longitud de los arreglos
    # será diferente, por lo que se agregará un 0 al final del segundo
    # para igualarlos
    if len(temp1) != len(temp2):
        temp2.append(0)
        
    # Se agregan los dos arreglos en una tabla
    tabla = []
    tabla.append(temp1)
    tabla.append(temp2)
    
    # Se itera atraves de las filas de la tabla y se calculan
    for i in range(2, len(polinomio)):
        tabla.append(fila(tabla, i))
        
    # Se agrega al inicio de cada fila el grado del polinomio asociado
    n = len(tabla)
    for i in range(n):
        tabla[i] = ["s**" + str(n - i - 1)] + tabla[i]
        
    # Se imprimen cada una de las lineas de la tabla
    if linea_cero:
        print "Se saco la derivada de una linea para sacar la siguiente,"
        print "debido a que es linealmente dependiente de la anterior."
        
    print tabulate(tabla)
```


```
A = [1, 2, 1, 2]
for r in roots(A):
    print r
routh_hurwitz(A)
```

    (-2+0j)
    (-1.38777878078e-16+1j)
    (-1.38777878078e-16-1j)
    ----  -  -
    s**3  1  1
    s**2  2  2
    s**1  2  0
    s**0  2  0
    ----  -  -



```
A = [1, 0, -3, 2]
for r in roots(A):
    print r
routh_hurwitz(A)
```


```
A = [1, 2, 24, 48, -25, -50]
for r in roots(A):
    print r
routh_hurwitz(A)
```

    (-1.11022302463e-16+5j)
    (-1.11022302463e-16-5j)
    (1+0j)
    (-2+0j)
    (-1+0j)
    ----  ---  ---  ---
    s**5    1   24  -25
    s**4    2   48  -50
    s**3    6   48    0
    s**2   32  -50    0
    s**1   57    0    0
    s**0  -50    0    0
    ----  ---  ---  ---



```
A = [1, 2, 3, 4, 5]
for r in roots(A):
    print r
routh_hurwitz(A)
```

    (0.287815479558+1.41609308017j)
    (0.287815479558-1.41609308017j)
    (-1.28781547956+0.857896758328j)
    (-1.28781547956-0.857896758328j)
    ----  --  -  -
    s**4   1  3  5
    s**3   2  4  0
    s**2   1  5  0
    s**1  -6  0  0
    s**0   5  0  0
    ----  --  -  -



```
A = [1, 3, 1, 5]
for r in roots(A):
    print r
routh_hurwitz(A)
```

    (-3.17998107216+0j)
    (0.0899905360791+1.25069504925j)
    (0.0899905360791-1.25069504925j)
    ----  --  -
    s**3   1  1
    s**2   3  5
    s**1  -1  0
    s**0   5  0
    ----  --  -



```

```

627ffaf940490c13130ec20a74b8b7a18c01f795
ipynb
Jupyter Notebook
IPythonNotebooks/Propedeutico/Clase 02.ipynb
chelizalde/DCA
34fd4d500117a9c0a75b979b8b0f121c1992b9dc
MIT
IPythonNotebooks/Propedeutico/Clase 02.ipynb
chelizalde/DCA
34fd4d500117a9c0a75b979b8b0f121c1992b9dc
MIT
IPythonNotebooks/Propedeutico/Clase 02.ipynb
chelizalde/DCA
34fd4d500117a9c0a75b979b8b0f121c1992b9dc
MIT
2021-03-20T12:44:13.000Z
2021-03-20T12:44:13.000Z
Qwen/Qwen-72B
1. YES
2. YES