2017-09-22 81 views
5

Próbuję obliczyć wartości własne symbolicznej złożonej macierzy M o rozmiarze 3x3. W niektórych przypadkach eigenvals() działa idealnie. Na przykład, następujący kod:Obliczanie symbolicznych wartości własnych za pomocą sympy

import sympy as sp 

kx = sp.symbols('kx') 
x = 0. 

M = sp.Matrix([[0., 0., 0.], [0., 0., 0.], [0., 0., 0.]]) 
M[0, 0] = 1. 
M[0, 1] = 2./3. 
M[0, 2] = 2./3. 
M[1, 0] = sp.exp(1j*kx) * 1./6. + x 
M[1, 1] = sp.exp(1j*kx) * 2./3. 
M[1, 2] = sp.exp(1j*kx) * -1./3. 
M[2, 0] = sp.exp(-1j*kx) * 1./6. 
M[2, 1] = sp.exp(-1j*kx) * -1./3. 
M[2, 2] = sp.exp(-1j*kx) * 2./3. 

dict_eig = M.eigenvals() 

zwraca mi 3 poprawnych złożone symboliczne wartości własne M. Jednakże kiedy ustawić x=1., pojawia się następujący błąd:

raise MatrixError("Could not compute eigenvalues for {}".format(self))

Próbowałem też do obliczania wartości własnych w następujący sposób:

lam = sp.symbols('lambda') 
cp = sp.det(M - lam * sp.eye(3)) 
eigs = sp.solveset(cp, lam) 

ale zwraca mi ConditionSet W każdym przypadku, nawet gdy eigenvals() puszka wykonuj pracę.

Czy ktoś wie, jak prawidłowo rozwiązać ten problem wartości własnej, dla dowolnej wartości x?

Odpowiedz

2

Twoja definicja M spowodowała, że ​​życie stało się zbyt trudne dla SymPy, ponieważ wprowadziło liczby zmiennoprzecinkowe. Kiedy chcesz mieć symboliczne rozwiązanie, unikaj pływaków. Oznacza to, że:

  • zamiast 1./3. (liczba zmiennoprzecinkowa Pythona) używać sp.Rational(1, 3) (numer racjonalnego SymPy), albo sp.S(1)/3, który ma ten sam efekt, ale jest łatwiejsze do pisania.
  • zamiast 1j (jednostka urojona Pythona) używać sp.I (jednostka urojona SymPy za)
  • zamiast x = 1., pisać x = 1 (Python 2.7 nawyków i SymPy pójść źle razem).

Z tych zmian albo solveset lub solve znaleźć wartości własne, choć solve dostaje je znacznie szybciej. Ponadto, można dokonać Poly obiekt i zastosować roots do tego, co jest prawdopodobnie najbardziej skuteczny: (. Byłoby łatwiej zrobić from sympy import * niż wpisać wszystkie te sp)

M = sp.Matrix([ 
    [ 
     1, 
     sp.Rational(2, 3), 
     sp.Rational(2, 3), 
    ], 
    [ 
     sp.exp(sp.I*kx) * sp.Rational(1, 6) + x, 
     sp.exp(sp.I*kx) * sp.Rational(1, 6), 
     sp.exp(sp.I*kx) * sp.Rational(-1, 3), 
    ], 
    [ 
     sp.exp(-sp.I*kx) * sp.Rational(1, 6), 
     sp.exp(-sp.I*kx) * sp.Rational(-1, 3), 
     sp.exp(-sp.I*kx) * sp.Rational(2, 3), 
    ] 
]) 
lam = sp.symbols('lambda') 
cp = sp.det(M - lam * sp.eye(3)) 
eigs = sp.roots(sp.Poly(cp, lam)) 


I "nie jestem całkiem pewien, dlaczego metoda EMSS Sympy zgłasza awarię, nawet pomimo powyższych modyfikacji. Jak widać in the source, nie robi to nic więcej niż powyższy kod: wywołaj roots na charakterystycznym wielomianie. Różnica wydaje się być w sposób wielomian ten jest tworzony: M.charpoly(lam) powraca

PurePoly(lambda**3 + (I*sin(kx)/2 - 5*cos(kx)/6 - 1)*lambda**2 + (-I*sin(kx)/2 + 11*cos(kx)/18 - 2/3)*lambda + 1/6 + 2*exp(-I*kx)/3, lambda, domain='EX') 

z tajemniczym (do mnie) domain='EX'. Następnie aplikacja roots zwraca {}, nie znaleziono żadnych korzeni. Wygląda na niedociągnięcie implementacji.

+0

Bardzo dziękuję za pomoc. Wydaje się, że mój problem pojawił się raczej w wyniku użycia 1j zamiast sp.I, ale użycie Rational z pewnością pomogło! Problem rozwiązany dla mnie, ale wciąż jest problem z eigenvalami SymPy ... – Azlof

+0

Uprościliśmy Twój przykład i opublikowałem [jako problem z SymPy] (https://github.com/sympy/sympy/issues/13340) – FTP

+0

Problem rozwiązany na github.Dla tych, którzy są zainteresowani, poprawka została pchnięta w gałęzi master programu SymPy. Dzięki Michelle! – Azlof