2016-03-29 27 views
4

Ben sympy ile polinom fonksiyonu çözerek sorunları var sympy. Aşağıdaki örnek, yönetemediğim bir hata iletisi veren bir durumu gösterir. Polinom basitleşirse çözücü düzgün çalışır. Lütfen sisteminizdeki hatayı kontrol etmek için kodu kopyalayıp yapıştırın.Polinom fonksiyon

import sympy 
from sympy import I 
omega = sympy.symbols('omega') 

def function(omega): 
    return - 0.34*omega**4  \ 
      + 7.44*omega**3  \ 
      + 4.51*I*omega**3 \ 
      + 87705.64*omega**2 \ 
      - 53.08*I*omega**2 \ 
      - 144140.03*omega \ 
      - 22959.95*I*omega \ 
      + 42357.18 + 50317.77*I 
sympy.solve(function(omega), omega) 

Bir sonuca nasıl ulaşabileceğimi biliyor musunuz? Yardımlarınız için teşekkür ederiz.

DÜZENLEME:

Bu hata iletisi:


TypeError         Traceback (most recent call last) 
<ipython-input-7-512446a62fa9> in <module>() 
     1 def function(omega): 
     2  return - 0.34*omega**4     + 7.44*omega**3     + 4.51*I*omega**3    + 87705.64*omega**2    - 53.08*I*omega**2    - 144140.03*omega    - 22959.95*I*omega    + 42357.18 + 50317.77*I 
----> 3 sympy.solve(function(omega), omega) 

C:\Anaconda\lib\site-packages\sympy\solvers\solvers.pyc in solve(f, *symbols, **flags) 
    1123  # restore floats 
    1124  if floats and solution and flags.get('rational', None) is None: 
-> 1125   solution = nfloat(solution, exponent=False) 
    1126 
    1127  if check and solution: # assumption checking 

C:\Anaconda\lib\site-packages\sympy\core\function.pyc in nfloat(expr, n, exponent) 
    2463    return type(expr)([(k, nfloat(v, n, exponent)) for k, v in 
    2464        list(expr.items())]) 
-> 2465   return type(expr)([nfloat(a, n, exponent) for a in expr]) 
    2466  rv = sympify(expr) 
    2467 

C:\Anaconda\lib\site-packages\sympy\core\function.pyc in nfloat(expr, n, exponent) 
    2497  return rv.xreplace(Transform(
    2498   lambda x: x.func(*nfloat(x.args, n, exponent)), 
-> 2499   lambda x: isinstance(x, Function))) 
    2500 
    2501 

C:\Anaconda\lib\site-packages\sympy\core\basic.pyc in xreplace(self, rule) 
    1085 
    1086   """ 
-> 1087   value, _ = self._xreplace(rule) 
    1088   return value 
    1089 

C:\Anaconda\lib\site-packages\sympy\core\basic.pyc in _xreplace(self, rule) 
    1093   """ 
    1094   if self in rule: 
-> 1095    return rule[self], True 
    1096   elif rule: 
    1097    args = [] 

C:\Anaconda\lib\site-packages\sympy\core\rules.pyc in __getitem__(self, key) 
    57  def __getitem__(self, key): 
    58   if self._filter(key): 
---> 59    return self._transform(key) 
    60   else: 
    61    raise KeyError(key) 

C:\Anaconda\lib\site-packages\sympy\core\function.pyc in <lambda>(x) 
    2496 
    2497  return rv.xreplace(Transform(
-> 2498   lambda x: x.func(*nfloat(x.args, n, exponent)), 
    2499   lambda x: isinstance(x, Function))) 
    2500 

C:\Anaconda\lib\site-packages\sympy\core\function.pyc in nfloat(expr, n, exponent) 
    2463    return type(expr)([(k, nfloat(v, n, exponent)) for k, v in 
    2464        list(expr.items())]) 
-> 2465   return type(expr)([nfloat(a, n, exponent) for a in expr]) 
    2466  rv = sympify(expr) 
    2467 

C:\Anaconda\lib\site-packages\sympy\core\function.pyc in nfloat(expr, n, exponent) 
    2463    return type(expr)([(k, nfloat(v, n, exponent)) for k, v in 
    2464        list(expr.items())]) 
-> 2465   return type(expr)([nfloat(a, n, exponent) for a in expr]) 
    2466  rv = sympify(expr) 
    2467 

TypeError: __new__() takes exactly 3 arguments (2 given) 
+0

ve bu hata mesajı nedir "Ben yönetemez bir hata mesajı verir"? –

+0

SO'ya hoş geldiniz! Aynı soruna sahip olan kullanıcıların bunu google'a girip, umarız çözülecek sorunuzu bulabilmeleri için hata mesajını sorunuza eklemelisiniz. –

+0

Hata mesajının son satırını ekledim. Umarım bu yardımcı olur. – Aloex

cevap

6

ben sympy aşina değilim yorumlarda belirtildiği gibi, ama burada denklemin köklerini bulmak için nasıl keyfi hassasiyet mpmath modülünü kullanarak. Python yüzen ile hassas sorunları önlemek amacıyla

, onu tamsayılar onları oluşturmak için elverişli olmadığı sürece, dize şeklinde mpmath için yüzen geçmek zamanki bu.

from mpmath import mp 

I = mp.mpc(0, 1) 

def func(x): 
    return (-mp.mpf('0.34') * x ** 4 
     + mp.mpf('7.44') * x ** 3 
     + mp.mpf('4.51') * I * x ** 3 
     + mp.mpf('87705.64') * x ** 2 
     - mp.mpf('53.08') * I * x ** 2 
     - mp.mpf('144140.03') * x 
     - mp.mpf('22959.95') * I * x 
     + mp.mpf('42357.18') + mp.mpf('50317.77') * I) 

mpf keyfi kesinlikte: Ben burada

mpmath sözdizimi doğrudan tercüme senin denklem ... senin katsayıları, oldukça düşük hassasiyet var ama yine de çünkü, gerçekten Denkleminizi ile ilgili bir sorun değil sanırım kurucu şamandıra, mpc karmaşık sayı kurucudur. Bu kurucuları nasıl arayacağınıza dair bilgi için lütfen mpmath belgelerine bakın.

Ancak biz böyle I hakkında karışıklık gerekmez: sadece karmaşık sayılar olarak doğrudan katsayıları tanımlayabilirsiniz.

from __future__ import print_function 
from mpmath import mp 

# set precision to ~30 decimal places 
mp.dps = 30 

def func(x): 
    return (mp.mpf('-0.34') * x ** 4 
     + mp.mpc('7.44', '4.51') * x ** 3 
     + mp.mpc('87705.64', '-53.08') * x ** 2 
     + mp.mpc('-144140.03', '-22959.95') * x 
     + mp.mpc('42357.18', '50317.77')) 

x = mp.findroot(func, 1) 
print(x) 
print('test:', func(x)) 

çıkış

(1.35717989161653180236360985534 - 0.202974596285109153971324419197j) 
test: (-3.2311742677852643549664402034e-26 + 6.4623485355705287099328804068e-27j) 

Ama nasıl diğer kökleri bulabilirim? Basit!

u f (x) bir kök olmak edelim. Sonra f (x) = g (x) (x - u) ve g (x)'un herhangi bir kökü de f (x)'un köküdür. Biz uygun başka listede bu yeni işlevi depolanması, bir listeye her bulundu kökünü kaydeder ve daha sonra önceki işlevinden yeni bir fonksiyon inşa eden for döngü kullanarak bu birden çok kez yapabilirsiniz. önerilir gibi karmaşık kökleri ararken Bu versiyonda

, ben, "Muller" çözücüsü kullanmak, ama aslında çözücüsü "sekant" varsayılan kullanarak aynı cevapları veriyor.

from __future__ import print_function 
from mpmath import mp 

# set precision to ~30 decimal places 
mp.dps = 30 

def func(x): 
    return (mp.mpf('-0.34') * x ** 4 
     + mp.mpc('7.44', '4.51') * x ** 3 
     + mp.mpc('87705.64', '-53.08') * x ** 2 
     + mp.mpc('-144140.03', '-22959.95') * x 
     + mp.mpc('42357.18', '50317.77')) 

x = mp.findroot(func, 1) 
print(x) 
print('test:', func(x)) 

funcs = [func] 
roots = [] 

#Find all 4 roots 
for i in range(4): 
    x = mp.findroot(funcs[i], 1, solver="muller") 
    print('%d: %s' % (i, x)) 
    print('test: %s\n%s\n' % (funcs[i](x), funcs[0](x))) 
    roots.append(x) 
    funcs.append(lambda x,i=i: funcs[i](x)/(x - roots[i]))  

lambda x,i=i: funcs[i](x)/(x - roots[i]) 

yılında çıktı

(1.35717989161653180236360985534 - 0.202974596285109153971324419197j) 
test: (-3.2311742677852643549664402034e-26 + 6.4623485355705287099328804068e-27j) 
0: (1.35717989161653180236360985534 - 0.202974596285109153971324419197j) 
test: (-3.2311742677852643549664402034e-26 + 6.4623485355705287099328804068e-27j) 
(-3.2311742677852643549664402034e-26 + 6.4623485355705287099328804068e-27j) 

1: (0.2859569222439674364374376897 + 0.465618100581394597267702975072j) 
test: (2.70967991111831205485691272044e-27 - 4.34146435347156317045282996313e-27j) 
(0.0 + 6.4623485355705287099328804068e-27j) 

2: (-497.86487129703641182172086688 + 6.49836193448774263077718499855j) 
test: (1.25428695883356196194609388771e-26 - 3.46609896266051486795576778151e-28j) 
(3.11655159180984988723836070362e-21 - 1.65830325771275337225587644119e-22j) 

3: (518.104087424352383171155113452 + 6.50370044356891310239702468087j) 
test: (-4.82755073209873082528016484574e-30 + 7.95804877623117526215e-31j) 
(-1.31713649147437845007238988587e-21 + 1.68350641700147843422461467477e-22j) 

bir varsayılan anahtar kelime argüman olarak i belirtmek böylece tanımlandı işlevi kullanıldığında i olduğunu değer. Aksi halde, geçerli değeri i, işlev çağrıldığında kullanılır; bu bizim istediğimiz değildir.


Birden çok kök bulma yöntemi, rastgele işlevler için kullanılabilir. Ancak, polinomları çözmek istediğimizde, mpmath'in tüm kökleri aynı anda bulabilmesi için daha iyi bir yolu vardır: polyroots işlevi. Sadece polinom katsayılarının bir listesini (veya bir kopyasını) vermeliyiz. Gördüğünüz gibi

from __future__ import print_function 
from mpmath import mp 

# set precision to ~30 decimal places 
mp.dps = 30 

coeff = (
    mp.mpf('-0.34'), 
    mp.mpc('7.44', '4.51'), 
    mp.mpc('87705.64', '-53.08'), 
    mp.mpc('-144140.03', '-22959.95'), 
    mp.mpc('42357.18', '50317.77'), 
) 

roots = mp.polyroots(coeff) 
for i, x in enumerate(roots): 
    print('%d: %s' % (i, x)) 
    print('test: %s\n' % mp.polyval(coeff, x)) 

çıkış

0: (1.35717989161653180236360985534 - 0.202974596285109153971324419197j) 
test: (6.4623485355705287099328804068e-27 - 6.4623485355705287099328804068e-27j) 

1: (0.2859569222439674364374376897 + 0.465618100581394597267702975072j) 
test: (0.0 + 0.0j) 

2: (-497.86487129703641182172086688 + 6.49836193448774263077718499855j) 
test: (-2.27689218423463552142807161949e-21 + 7.09242751778865525915133624646e-23j) 

3: (518.104087424352383171155113452 + 6.50370044356891310239702468087j) 
test: (-7.83663157514495734538720675411e-22 - 1.08373584941517766465574404422e-23j) 

, sonuçları önceki teknik ile elde edilenlere çok yakındır. polyroots kullanmak sadece daha doğru değil, kök için bir başlangıç ​​yaklaşımı belirlememize gerek duymamanın büyük avantajına sahiptir, mpmath kendisi için bir tane oluşturmak için yeterince akıllıdır.

+0

Çok teşekkürler! ... sry, diğer üç kökten bahsetmeyi özledim. :) Son olarak, kodunuzu senaryona koydum ve güzel bir çözüm buldum. Başlangıçta bir semp fonksiyonu ile gelirseniz tam bir örnek göndereceğim. – Aloex

+0

@Aloex: Semptomu kullanmasa bile cevabımın size yardımcı olduğuna sevindim. Lütfen [kabul etmeyi] düşünün (http://meta.stackexchange.com/questions/5234/how-does-accepting-an-answer-work). –

+0

Ayrıca 'mpmath.findroot' etrafında bir paket olan' sympy.nsolve' kullanabilirsiniz. – asmeurer

1

PM2Ring'in yardımıyla, dördüncü dereceden bir sıralama polinomunun verilişindeki katsayıları çıkartan ve çözen tam bir kod örneği oluşturdum.

import sympy 
from sympy import I 
import mpmath as mp 
import numpy as np 

omega = sympy.symbols('omega') 

# Define polynomial in sympy 
def function(omega): 
    return (- 0.34     *omega**4  \ 
      + (7.44 + 4.51*I)   *omega**3  \ 
      + (87705.64 - 53.08*I) *omega**2 \ 
      - (144140.03 + 22959.95*I)*omega \ 
      + 42357.18 + 50317.77*I) 

# Show defined polynomial 
print''+ str(function(omega).expand()) +'\n' 
result = sympy.Poly(function(omega).expand(), omega) 

# Obtain coefficients of the polynomial 
coeff = (
    mp.mpc(str(np.real(sympy.lambdify((I),result.coeffs()[0])(1j))) , str(np.imag(sympy.lambdify((I),result.coeffs()[0])(1j)))), 
    mp.mpc(str(np.real(sympy.lambdify((I),result.coeffs()[1])(1j))) , str(np.imag(sympy.lambdify((I),result.coeffs()[1])(1j)))), 
    mp.mpc(str(np.real(sympy.lambdify((I),result.coeffs()[2])(1j))) , str(np.imag(sympy.lambdify((I),result.coeffs()[2])(1j)))), 
    mp.mpc(str(np.real(sympy.lambdify((I),result.coeffs()[3])(1j))) , str(np.imag(sympy.lambdify((I),result.coeffs()[3])(1j)))), 
    mp.mpc(str(np.real(sympy.lambdify((I),result.coeffs()[4])(1j))) , str(np.imag(sympy.lambdify((I),result.coeffs()[4])(1j)))), 
) 

# Calculate roots of the polynomial 
roots = mp.polyroots(coeff) 
for no, frequency in enumerate(roots): 
    print frequency 

Çıktı

-0.34*omega**4 + 7.44*omega**3 + 4.51*I*omega**3 + 87705.64*omega**2 - 53.08*I*omega**2 - 144140.03*omega - 22959.95*I*omega + 42357.18 + 50317.77*I 

(1.35717989161653 - 0.202974596285109j) 
(0.285956922243967 + 0.465618100581395j) 
(-497.864871297036 + 6.49836193448774j) 
(518.104087424352 + 6.50370044356891j)