2015-02-06 21 views
5

Ho un'equazione polinomiale di 4 ° ordine e ho bisogno di trovare tutte le radici. Esempio semplice:SymPy non può risolvere l'equazione polinomiale del 4 ° ordine

from sympy import (Symbol,solve,I) 

a=4+5*I; b=3+7*I; c=12-56*I; d=33+56*I; e=345-67*I; x=Symbol('x') 
eq=a*x**4 + b*x**3 + c*x**2 + d*x +e 
solve(eq,x) 

Se a, b, c, d, e sono pure vero, allora funzionare bene. Ma nel mio caso tutti sono numeri complessi. Poi ho avuto chiamare:

PolynomialError: 'cannot return general quartic solution' 

trovo tipo di problema simile, e implementare la correzione: Description of the issue. Fix of the issue

ma in realtà non aiuta. V'è una sorta di strano problema, come ora la chiamata è (come modificato nella correzione):

PolynomialError: Cannot determine if `-((12 - 56*I)/(4 + 5*I) - 3*(3 + 7*I)**2/(8*(4 + 5*I)**2))**2/12 + (3 + 7*I)*((33 + 56*I)/(4*(4 + 5*I)) + (3 + 7*I)*(3*(3 + 7*I)**2/(256*(4 + 5*I)**2) - (12 - 56*I)/(16*(4 + 5*I)))/(4 + 5*I))/(4 + 5*I) - (345 - 67*I)/(4 + 5*I)` is nonzero. 

Ma per determinare se espressione di cui sopra è diverso da zero è la cosa più semplice, quindi non so dove il problema potrebbe essere.

risposta

1

Aggiornamento alla versione più recente di SymPy che supporta soluzioni quartiche arbitrarie.

+0

Dopo essere passato a sympy 0.7.6 (usando Python 3.4) v'è diverso massaggio errore : riga 103 in __nonzero__: genera TypeError ("non può determinare il valore di verità di \ n% s"% self) – K4stan

+0

Sembra che questo sia un bug, che è stato corretto nella versione git di SymPy. – asmeurer

+0

Sto usando la distribuzione anaconda di python. È possibile aggiornare in qualche modo a questa versione git lì, o devo riscrivere alcuni script? – K4stan

0

Se si desidera una soluzione più flessibile, si potrebbe risolvere per x usando la ricerca binaria con qualcosa come il seguente:

def maybeRightX(maybeX, polys): 
    sum = 0 
    for i in range(len(polys)): 
     sum += polys[i]*(maybeX ** i) 
    return sum 

def solve(y, polys): 
    lo = 0 
    hi = y 
    while lo <= hi: 
     mid = (lo + hi)//2 
     if (maybeRightX(mid, polys)) < y: 
      lo = mid + 1 
     else: 
      hi = mid - 1 
    return (hi + 1)