In [1]:

    

import numpy as np
import sympy as sy


    

In [2]:

    

s = sy.symbols('s', real=False)


    

In [3]:

    

np.roots?


    

In [5]:

    

poles=np.roots([1, 25.9, 168])


    

In [6]:

    

np.abs(poles)


    

    
Out[6]:





array([ 12.9614814,  12.9614814])

In [7]:

    

np.angle(poles)


    

    
Out[7]:





array([ 3.09949895, -3.09949895])

In [8]:

    

poles


    

    
Out[8]:





array([-12.95+0.54543561j, -12.95-0.54543561j])

In [9]:

    

G = 168/(s*(s**2 + 25.9*s + 168))


    

In [10]:

    

Y = sy.apart(G/s)
sy.pretty_print(Y)


    

    




1.0⋅(0.154166666666667⋅s + 2.99291666666667)   - -0.154166666666667    1.0
──────────────────────────────────────────── - ───────────────────── + ───
               2                                         s               2
          1.0⋅s  + 25.9⋅s + 168.0                                       s 

In [11]:

    

Y


    

    
Out[11]:





1.0*(0.154166666666667*s + 2.99291666666667)/(1.0*s**2 + 25.9*s + 168.0) - 0.154166666666667/s + 1.0/s**2

In [22]:

    

a = 12.95
omega = 0.54543561
b1 = 0.15416666667
rest = (2.9929166666667/b1-a)
b2=rest*b1/omega
print b1
print b2


    

    




0.15416666667
1.82690369866

In [23]:

    

den = s**2 + 25.9*s + 168
sy.together(b1*(s+a)/den + b2*omega/den)


    

    
Out[23]:





(0.15416666667*s + 2.9929166666667)/(s**2 + 25.9*s + 168)

In [32]:

    

zeta = 12.95/np.sqrt(0.5454**2 + 12.95**2)
phi = np.arccos(zeta)
np.exp(phi/np.tan(phi))


    

    
Out[32]:





2.7166768328237918

In [33]:

    

print zeta


    

    




0.999114306497

In [34]:

    

print phi


    

    




0.0420909556946

In [35]:

    

2.72/4


    

    
Out[35]:





0.68

In [38]:

    

omega0=np.sqrt(0.5454**2 + 12.95**2)


    

In [39]:

    

2.72/4/omega0


    

    
Out[39]:





0.052463145051601344

In [40]:

    

2.72/10/omega0


    

    
Out[40]:





0.020985258020640538

In [42]:

    

h=0.025
dreal = np.exp(h*(-12.95))*np.cos(0.5454*h)
dim = np.exp(h*(-12.95))*np.sin(0.5454*h)
print dreal
print dim
dpole = dreal + 1j*dim


    

    




0.723363830883
0.00986367710434

In [43]:

    

print np.abs(dpole)
print np.angle(dpole)


    

    




0.723431077544
0.013635

In [ ]: