使用 SymPy 查找复杂多项式的所有根答案

【问题标题】：Finding all roots of a complex polynomial with SymPy使用 SymPy 查找复杂多项式的所有根
【发布时间】：2013-04-05 09:31:19
【问题描述】：

我正在尝试用 SymPy 象征性地求解具有复数及其共轭的多项式。我想我已经走了很长一段路，但solve 没有给我任何解决方案，尽管多项式是可解的。

from sympy import *

# set up symbols
a, b = symbols("a b", real=True)
t = a+I*b
T = functions.conjugate(t)

# set up polynomial
a1=0.005+I*0.0009
a2=0.9+I*-0.9
a3=0.4+I*0.5
a4=8+I*-80
a5=284+I*-1.5
a6=27100+I*-11500
poly=t**2 * T * a1 + t * T * a2 + t**2 * a3 + T * a4 + t * a5 + a6

# Trying to solve symbolically...
solve([re(poly), im(poly)], a, b)
#     Output: []

# Solving numerically works, but only finds one solution...
nsolve((re(poly), im(poly)), (a, b), (0, 0))
#     Output: matrix(
#             [['-137.962596090596'],
#              ['52.6296963395752']])

# verify with two solutions obtained in Maxima
poly.subs({a:-137.9625935162095, b:52.6296992481203}).n()
#     Output: 0.000540354631040322 + 0.00054727003909351*I
poly.subs({a:-332.6474382554614+I*-185.9848818313149, b:258.0065640091016+I*-272.3344263478699}).n()
#     Output: -6.55448222470652e-12 - 1.41238056784605e-12*I

有什么想法吗？

【问题讨论】：

这确实很奇怪。只是一个需要确定的问题：鉴于第一个替代后的值远非零，是否可能无法解决 poly 并且最大值为您提供了数值近似？（我不太了解最大值，所以我不太清楚它的行为）
@EnricoGiampieri 很可能是 Maxima 给了我一个数字近似值，但在它的输出中并没有提到它。我也不是 Maxima 专家 :-( 无论如何，如果 SymPy 能给出 Maxima 给出的解决方案就好了，即使它们只是数字近似值。
我尝试了一些处理这种方程，看起来问题出在 sympy 求解机制中，它不能很好地处理复数。我试图找出我们是否错过了使其工作所需的某些条件......
这在我看来不像一个多项式。多项式具有a_0 z^0 + a_1 z^1 + ... 类型的常数系数，您所拥有的是a_1 z^2 z* 之类的项，即a1 项乘以它的复数conj。如果你有一个来自 sympy 的多项式 mpmath 应该会很好用 docs.sympy.org/dev/modules/mpmath/calculus/polynomials.html
SymPy 似乎无法解决这个问题，但即使它可以，它也不会为您提供 Maxima 的结果，因为您将 a 和 b 设置为真实的，所以solve() 只会给你a 和b 的真实值。

标签： python math polynomial-math sympy

【解决方案1】：

一个问题可能是您的系数包含浮点数。这通常不适用于符号软件。

使用 f=10000*simplify(re(poly)) 和 g=10000*simplify(im(poly)) 并编辑结果可得到具有整数系数的多项式。然后，CAS（在我的例子中是 Magma）可以生成 f 和 g 理想的三角形表示，它以多项式形式给出

a - 2483798807340123709959247/13545514719183259347425004828125*b^4
  + 66732206412048596386372283/541820588767330373897000193125*b^3 
  - 3849759933277117021785191063/86691294202772859823520030900*b^2
  + 9245906471290310401430681453/1733825884055457196470400618*b
  - 31414499425567273751868164900/866912942027728598235200309,

b^5 - 189465979625/206648369*b^4 
    + 330827538698125/826593476*b^3
    - 17645868534640625/206648369*b^2 
    + 1724106750659765625/206648369*b
    - 52548859891484375000/206648369

这告诉我们存在 5 个解决方案。第二个多项式的数值解是

174.10461010682254746847015187264557067610513554291323564643772 
 + 63.402741884833821878468926640811609033267039765104756747285816*i

174.104610106822547468470151872645570676105135542913235646437738
 - 63.402741884833821878468926640811609033267039765104756747285804*i

258.006564009101655109715962546854008929462784282347754971379392
  + 272.334426347869856080204881056671278679761260094680345276069337*i

258.006564009101655109715962546854008929462784282347754971379382
 - 272.334426347869856080204881056671278679761260094680345276069359*i

52.62969633957523864698147679803879873180829265956656342643376

产生了一个真正的解决方案。 sympy的数值结果正确且完整。

【讨论】：