Revision history [back]

EDIT: Some more info; I'm studying a dynamical system at an equilibrium point. I would like to get an expression for the eigenvectors and values assosiated with the center eigenspace of a "hopf bifurcation". Here is the code in question:

var('alpha1, alpha2, beta1, beta2, mu1, mu2, phi1, phi2, c1, c2, k1, k2, t ') ;
# suggested parameter values:

alpha1 = 0.7; alpha2 = 0.3
beta1 = 0.003; beta2 = 0.0015
mu1 = 0.6; mu2 = 0.4
phi1 = 2.3; phi2 = 4
c1 = 400000; c2 = 100
k1 = 0.0033100000000000;  
k2 = 0.000988067680992286;

# The vector field and its jacobian matrix:
f(P, M, L, G,k2) = [phi1*M*(1-M/c1) - (alpha1 + beta1)* P  - k1*P*G, alpha1*P - mu1*M, phi2*G*(1-G/c2) - (alpha2 + beta2)* L  + k2*P*G, alpha2*L - mu2*G];
J = jacobian(f, (P, M, L, G));

# Equilibrium solution:
M4 =((alpha1*alpha2 + alpha1*beta2)*c1*c2*k1*mu1*mu2 + (alpha1^2*alpha2*c1*phi1 - (alpha1*alpha2*c1*c2*k1 + (alpha1^2*alpha2 + alpha1*alpha2*beta1)*c1)*mu1)*phi2)/(alpha2*c1*c2*k1*k2*mu1^2 + alpha1^2*alpha2*phi1*phi2);

P4 = ((alpha2 + beta2)*c1*c2*k1*mu1^2*mu2 + (alpha1*alpha2*c1*mu1*phi1 - (alpha2*c1*c2*k1 + (alpha1*alpha2 + alpha2*beta1)*c1)*mu1^2)*phi2)/(alpha2*c1*c2*k1*k2*mu1^2 + alpha1^2*alpha2*phi1*phi2);

L4 = -((alpha1*alpha2 + alpha2*beta1)*c1*c2*k2*mu1^2*mu2 - alpha1^2*alpha2*c2*mu2*phi1*phi2 - (alpha1*alpha2*c1*c2*k2*mu1*mu2 - (alpha1^2*alpha2 +alpha1^2*beta2)*c2*mu2^2)*phi1)/(alpha2^2*c1*c2*k1*k2*mu1^2 + alpha1^2*alpha2^2*phi1*phi2);

G4 = -((alpha1*alpha2 + alpha2*beta1)*c1*c2*k2*mu1^2 - alpha1^2*alpha2*c2*phi1*phi2 -(alpha1*alpha2*c1*c2*k2*mu1 - (alpha1^2*alpha2 + alpha1^2*beta2)*c2*mu2)*phi1)/(alpha2*c1*c2*k1*k2*mu1^2 + alpha1^2*alpha2*phi1*phi2);

# Now I evaluate it at the fixed point:
J4 = J(P=P4, G=G4, M=M4, L=L4)
J4 = J4(k2=k2)

Now if I do something like;

J4.eigenvalues()

I get an lengthy radical expression. And when computing the eigenvectors:

D, P = J4.eigenmatrix_left()

I get an error...

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "_sage_input_129.py", line 10, in <module>
  exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("RCwgUCA9IEo0LmVpZ2VubWF0cml4X2xlZnQoKQ=="),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))
  File "", line 1, in <module>

  File "/tmp/tmpYDeAbA/___code___.py", line 2, in <module>
  exec compile(u'D, P = J4.eigenmatrix_left()
  File "", line 1, in <module>

  File "sage/matrix/matrix2.pyx", line 5617, in sage.matrix.matrix2.Matrix.eigenmatrix_left (build/cythonized/sage/matrix/matrix2.c:39965)
  IndexError: list index out of range