[Numpy-discussion] Perturbations in Polynomial roots

31 Jul 2019

      Dear all,

I am not sure if I misunderstood the Polynomial.roots() method, but I 
get strange behavior. The roots are way outside of the expected values. 
That is not the case with np.roots(). What am I doing wrong?

Code:

import numpy as np

def roots1(coef, x, y):
     for i in range(y.size):
         coefN = coef.copy()
         coefN[0] -= y[i]
         print(np.roots(coefN[::-1]) - x[i])

def roots2(coef, x, y):
     for i in range(y.size):
         print((np.polynomial.polynomial.Polynomial(coef) - y[i]).roots() 
- x[i])

x = np.arange(4) * 0.1
alpha, beta, gamma = 0.16, 1, 1e-15
coef = np.array([alpha, beta, gamma])
y = alpha + beta*x + gamma*x**2
print("roots")
roots1(coef, x, y)
print('Polynomial')
roots2(coef, x, y)

Outputs:

roots
[-1.e+15  0.e+00]
[-1.e+15  0.e+00]
[-1.00000000e+15  2.77555756e-17]
[-1.00000000e+15  5.55111512e-17]
Polynomial
[-1.e+15  0.e+00]
[-1.0e+15  2.5e-02]
[-1.0e+15 -7.5e-02]
[-1.e+15 -5.e-02]

Best Regards,
Cédric Hannotier