In [1]:

    

from sympy import *
init_session()


    

    




IPython console for SymPy 0.7.4.1-git (Python 3.3.3-64-bit) (ground types: python)

These commands were executed:
>>> from __future__ import division
>>> from sympy import *
>>> x, y, z, t = symbols('x y z t')
>>> k, m, n = symbols('k m n', integer=True)
>>> f, g, h = symbols('f g h', cls=Function)

Documentation can be found at http://www.sympy.org

In [4]:

    

var('a1:17 A M E')
B = S.One
f1 = lambda p: p**a1 / (a2**a1 + p**a1)
f2 = lambda p: a4 * a5**a3 / (a5**a3 + p**a3)
p = lambda E, a15: (a14 * E * B) / a15
dA = f1(p(E,a15))*(a6+a7*M)-a8*A-a9*A**2
dM = a10*f2(p(E,a15))*A-f1(p(E,a15))*a16*a7*M-a11*M
dE = a12*(1-f2(p(E,a15)))*A-a13*E


    

In [7]:

    

X = Matrix([dA, dM, dE])
Y = Matrix([A, M, E])
J = X.jacobian(Y)


    

In [8]:

    

J


    

    
Out[8]:





$$\left[\begin{matrix}- 2 A a_{9} - a_{8} & \frac{a_{7} \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}}{a_{2}^{a_{1}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}} & - \frac{a_{1} \left(\frac{E a_{14}}{a_{15}}\right)^{2 a_{1}} \left(M a_{7} + a_{6}\right)}{E \left(a_{2}^{a_{1}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}\right)^{2}} + \frac{a_{1} \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}} \left(M a_{7} + a_{6}\right)}{E \left(a_{2}^{a_{1}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}\right)}\\\frac{a_{10} a_{4} a_{5}^{a_{3}}}{a_{5}^{a_{3}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{3}}} & - a_{11} - \frac{a_{16} a_{7} \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}}{a_{2}^{a_{1}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}} & - \frac{A a_{10} a_{3} a_{4} a_{5}^{a_{3}} \left(\frac{E a_{14}}{a_{15}}\right)^{a_{3}}}{E \left(a_{5}^{a_{3}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{3}}\right)^{2}} + \frac{M a_{1} a_{16} a_{7} \left(\frac{E a_{14}}{a_{15}}\right)^{2 a_{1}}}{E \left(a_{2}^{a_{1}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}\right)^{2}} - \frac{M a_{1} a_{16} a_{7} \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}}{E \left(a_{2}^{a_{1}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{1}}\right)}\\a_{12} \left(- \frac{a_{4} a_{5}^{a_{3}}}{a_{5}^{a_{3}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{3}}} + 1\right) & 0 & \frac{A a_{12} a_{3} a_{4} a_{5}^{a_{3}} \left(\frac{E a_{14}}{a_{15}}\right)^{a_{3}}}{E \left(a_{5}^{a_{3}} + \left(\frac{E a_{14}}{a_{15}}\right)^{a_{3}}\right)^{2}} - a_{13}\end{matrix}\right]$$

In [16]:

    

substitutions = {
    a1: Integer(2),
    a2: Integer(2),
    a3: Integer(3),
    a4: Rational(7, 10),
    a5: Integer(1),
    a6: Rational(2, 100),
    a7: Integer(20),
    a8: Integer(1),
    a9: Integer(1),
    a10: Integer(1),
    a11: Rational(1, 100),
    a12: Rational(1, 10),
    a13: Rational(3, 10),
    a14: Integer(50),
    a16: Rational(1, 10)
}


    

In [18]:

    

Jprime = J.subs(substitutions)


    

In [21]:

    

[N(lambda_) for lambda_ in Jprime.subs({ A: 0.07016951, M: 0.007120757, E: 0.02117123, a15: 0.5708 }).eigenvals()]


    

    
Out[21]:





$$\begin{bmatrix}-2.2935285008276, & 5.77040297973451 \cdot 10^{-5} - 0.544873441891618 i, & 5.77040297973451 \cdot 10^{-5} + 0.544873441891618 i\end{bmatrix}$$

In [30]:

    

data = data[1:]
Jprime.subs({ A: data[0][1], M: data[0][2], E: data[0][3], a15: data[0][0] })


    

    
Out[30]:





$$\left[\begin{matrix}-1.037372 & 18.7593870418274 & 0.379702531800409\\0.0014849515486757 & -1.88593870418274 & -0.01390645344825\\0.0998515048451324 & 0 & -0.298664410993172\end{matrix}\right]$$

In [36]:

    

all_eigenvals = [[pp] + [N(lambda_) for lambda_ in Jprime.subs({ A: aa, M: mm, E: ee, a15: pp }).eigenvals()] for pp, aa, mm, ee in data]


    

In [37]:

    

all_eigenvals[:10]


    

    
Out[37]:





$$\begin{bmatrix}\begin{bmatrix}0.0399856, & -1.93593465232925 - 1.0 \cdot 10^{-23} i, & -1.01641960116789 + 3.0 \cdot 10^{-23} i, & -0.269620861678773 - 1.0 \cdot 10^{-23} i\end{bmatrix}, & \begin{bmatrix}0.0499764, & -0.262978566307137 - 3.0 \cdot 10^{-23} i, & -1.91968845355413 - 1.0 \cdot 10^{-23} i, & -0.965565865446056 + 3.0 \cdot 10^{-23} i\end{bmatrix}, & \begin{bmatrix}0.0599679, & -0.266425232119731 + 3.0 \cdot 10^{-21} i, & -0.873727911515919 + 1.0 \cdot 10^{-23} i, & -1.9191238296787 - 1.0 \cdot 10^{-23} i\end{bmatrix}, & \begin{bmatrix}0.0699633, & -0.730079208268977 + 1.0 \cdot 10^{-23} i, & -0.294215943602127 - 1.0 \cdot 10^{-23} i, & -1.93611806251275 + 2.0 \cdot 10^{-21} i\end{bmatrix}, & \begin{bmatrix}0.0799627, & -0.447411295713294 - 0.0586429669609455 i, & -0.447411295713294 + 0.0586429669609455 i, & -1.96563179691282\end{bmatrix}, & \begin{bmatrix}0.0899609, & -1.99964000258732, & -0.384432036741912 - 0.228390825577553 i, & -0.384432036741912 + 0.228390825577553 i\end{bmatrix}, & \begin{bmatrix}0.0999517, & -0.328868857637807 - 0.306154460708728 i, & -0.328868857637807 + 0.306154460708728 i, & -2.03209429833395\end{bmatrix}, & \begin{bmatrix}0.109932, & -2.06038127732021, & -0.282400540193982 + 0.356589298614659 i, & -0.282400540193982 - 0.356589298614659 i\end{bmatrix}, & \begin{bmatrix}0.119902, & -2.0841497652745, & -0.244333388562074 + 0.391774823394856 i, & -0.244333388562074 - 0.391774823394856 i\end{bmatrix}, & \begin{bmatrix}0.129863, & -2.1039565733311, & -0.213202589963214 - 0.417436578142293 i, & -0.213202589963214 + 0.417436578142293 i\end{bmatrix}\end{bmatrix}$$

In [38]:

    

import pickle


    

In [39]:

    

pickle.dump(data, open('data.pickle', 'wb'))
pickle.dump(all_eigenvals, open('eigenvals.pickle', 'wb'))


    

In [ ]: