为什么 numpy 的 root 失败了？

Question

为什么我的回答来自np.roots different from np.polynomial.polynomial.polyroots？

定义多项式

# a specific polynomial x**0 + x**1 + x**2 + x**3
p = [1, -2.8176255165067872, -0.97639120853458261, -0.86023870029448335]

这里有一个简单的例子来说明区别，

import numpy as np
r1 = np.roots(p); r2 = np.polynomial.polynomial.polyroots(p)
f = lambda x: np.sum([x**i*j for i,j in enumerate(p)])

print "{:>10} {:>10}".format("roots","polyroots")
for i,j in zip(r1, r2): # test should return 0
    print "{:10.5f} {:10.5f}".format(np.abs(f(i)),np.abs(f(j)))

输出明显不为零

     roots  polyroots
  46.41221    0.00000
   1.97595    0.00000
   1.97595    0.00000

正确的比较案例

相比之下，Mathematica 正确求根：

fn[x_] := 1.` - 2.817625516506788` x - 0.97639120853458261` x^2 - 0.8602387002944835` x^3
Roots[fn[x] == 0, x]

提供的根为：

x == -0.723475 - 1.78978 I || x == -0.723475 + 1.78978 I || x == 0.311926

测试验证了这一点：

fn[-0.7234748700272414` - 1.7897835665374093` I]
-4.44089*10^-16 - 2.66454*10^-15

Answer 1

numpy.polynomial 中的代码比 numpy.roots（和 numpy.poly1d 等）中的代码更新。在新的多项式代码中，系数顺序的约定发生了变化。在新代码中，系数是按递增顺序给出的，而在旧代码中，先给出最高阶的系数。

In [98]: p = [1, -2.8176255165067872, -0.97639120853458261, -0.86023870029448335]

In [99]: np.roots(p[::-1])
Out[99]: array([-0.72347487+1.78978357j, -0.72347487-1.78978357j,  0.31192616+0.j        ])

In [100]: np.polynomial.polynomial.polyroots(p)
Out[100]: array([-0.72347487-1.78978357j, -0.72347487+1.78978357j,  0.31192616+0.j        ])

为什么 numpy 的 root 失败了？

Why is numpy's root failing?

python

numpy

polynomial-math

polynomials

正确的比较案例