Корни полиномов Python неточны - PullRequest
1 голос
/ 01 ноября 2019

Алгоритм, который я пытаюсь реализовать, требует нахождения корней полинома 10-й степени, который я создал с помощью sympy, он выглядит следующим образом:

import sympy
import numpy as np
det = sympy.Poly(1.3339507303385e-16*z**10 + 6.75390067076469e-14*z**9 + 7.18791227134351e-12*z**8 + 2.27504286959352e-10*z**7 + 2.37058998324426e-8*z**6 + 1.63916629439745e-6*z**5 + 3.0608041245671e-5*z**4 + 4.83564348562906e-8*z**3 + 2.0248073519853e-5*z**2 - 4.73126177166988e-7*z + 1.1495883206077e-6)

Для нахождения корней полинома,Я использую следующий код:

coefflist = det.coeffs()
solutions = np.roots(coefflist)
print(coefflist)
[1.33395073033850e-16, 6.75390067076469e-14, 7.18791227134351e-12, 2.27504286959352e-10, 2.37058998324426e-8, 1.63916629439745e-6, 3.06080412456710e-5, 4.83564348562906e-8, 2.02480735198530e-5, -4.73126177166988e-7, 1.14958832060770e-6]

print(solutions)
[-3.70378229e+02+0.00000000e+00j -1.18366138e+02+0.00000000e+00j
  2.71097137e+01+5.77011644e+01j  2.71097137e+01-5.77011644e+01j
 -3.59084863e+01+1.44819591e-02j -3.59084863e+01-1.44819591e-02j
  2.60969082e-03+7.73805425e-01j  2.60969082e-03-7.73805425e-01j
  1.42936329e-02+2.49877948e-01j  1.42936329e-02-2.49877948e-01j]

Однако, когда я заменяю z корнем, скажем, первый, результат будет не нулевым, а некоторым числом:

print(det.subs(z,solutions[0]))
-1.80384169514123e-6

Я бы ожидал, что результат, вероятно, не целое число 0, но 1e-6 довольно плохо (должно быть ноль, верно?). Есть ли ошибка в моем коде? Это неточность нормально? Любые мысли / предложения будут полезны. Есть ли более точная альтернатива для вычисления корней полинома 10-й степени?

1 Ответ

0 голосов
/ 01 ноября 2019

Вам не нужно сочувствовать, методов в numpy вполне достаточно. Определите полином по его списку коэффициентов и вычислите корни

p=[1.33395073033850e-16, 6.75390067076469e-14, 7.18791227134351e-12, 2.27504286959352e-10, 2.37058998324426e-8, 1.63916629439745e-6, 3.06080412456710e-5, 4.83564348562906e-8, 2.02480735198530e-5, -4.73126177166988e-7, 1.14958832060770e-6]
sol= np.roots(p); sol

, дающие результат

array([ -3.70378229e+02 +0.00000000e+00j,  -1.18366138e+02 +0.00000000e+00j,
         2.71097137e+01 +5.77011644e+01j,   2.71097137e+01 -5.77011644e+01j,
        -3.59084863e+01 +1.44819592e-02j,  -3.59084863e+01 -1.44819592e-02j,
         2.60969082e-03 +7.73805425e-01j,   2.60969082e-03 -7.73805425e-01j,
         1.42936329e-02 +2.49877948e-01j,   1.42936329e-02 -2.49877948e-01j])

и оцените полиномы по этим приблизительным корням

np.polyval(p,sol)

присвоение массива

array([  2.28604877e-06 +0.00000000e+00j,   1.30435230e-10 +0.00000000e+00j,
         1.05461854e-11 -7.56043461e-12j,   1.05461854e-11 +7.56043461e-12j,
        -3.98439686e-14 +6.84489332e-17j,  -3.98439686e-14 -6.84489332e-17j,
         1.18584613e-20 +1.59976730e-21j,   1.18584613e-20 -1.59976730e-21j,
         6.35274710e-22 +1.74700545e-21j,   6.35274710e-22 -1.74700545e-21j])

Очевидно, что оценка полинома, близкого к корню, включает в себя множество катастрофических отмен, то есть промежуточные члены имеют большие значения противоположного знака и удаляются, но их ошибки пропорциональны их исходнымразмеры. Чтобы получить оценку объединенного размера ошибки, замените полиномиальные коэффициенты их абсолютными значениями, а также точками оценки.

np.polyval(np.abs(p),np.abs(sol))

, что приведет к

array([  1.81750423e+10,   8.40363409e+05,   
         8.08166359e+03,   8.08166359e+03,
         2.44160616e+02,   2.44160616e+02,
         2.50963696e-05,   2.50963696e-05,
         2.65889696e-06,   2.65889696e-06])

В случае первогоroot, масштаб, умноженный на машинную константу, дает масштаб ошибки 1e+10*1e-16=1e-6, что означает, что значение в корне равно нулю в рамках плавающей запятой двойной точности.

...