2016-10-29 9 views
0

Итак, я в конечном счете пытаюсь использовать правило Хорнера (http://mathworld.wolfram.com/HornersRule.html) для вычисления полиномов и создания функции для вычисления полинома. Как бы то ни было, моя проблема связана с тем, как я написал функцию; он работает для простых многочленов типа 3x^2 + 2x^1 + 5 и т. д. Но как только вы получите оценку полинома с числом с плавающей запятой (что-то сумасшедшее, например, 1.8953343e-20 и т. Д.), Оно теряет свою точность. Поскольку я использую эту функцию для оценки корней полинома с использованием метода Ньютона (http://tutorial.math.lamar.edu/Classes/CalcI/NewtonsMethod.aspx), мне нужно это, чтобы быть точным, поэтому оно не теряет своей ценности за счет небольшой ошибки округления и еще чего-то.Прецизионная ошибка в Python

У меня уже есть проблемы с двумя другими людьми, что проблема лежит в функции calculatePoly(), а не на моих других функциях, которые оценивают метод Ньютона. Кроме того, я первоначально оценивал полином обычно (умножая х на степень, умножая на константу и т. Д.), И он вытащил правильный ответ. Тем не менее, назначение требует использования правила Хорнера для упрощения вычислений.

Это мой следующий код:

def evaluatePoly(poly, x_): 
    """Evaluates the polynomial at x = x_ and returns the result as a floating 
point number using Horner's rule""" 
    #http://mathworld.wolfram.com/HornersRule.html 
    total = 0. 
    polyRev = poly[::-1] 
    for nn in polyRev: 
     total = total * x_ 
     total = total + nn 
    return total 

Примечание: Я уже попытался установить Nn, x_, (всего * Х-) в качестве поплавков с помощью поплавка().

Это выход я получаю:

Polynomial: 5040x^0 + 1602x^1 + 1127x^2 - 214x^3 - 75x^4 + 4x^5 + 1x^6 

Derivative: 1602x^0 + 2254x^1 - 642x^2 - 300x^3 + 20x^4 + 6x^5 

(6.9027369297630505, False) 

Starting at 100.00, no root found, last estimate was 6.90, giving value f(6.90) = -6.366463e-12 
(6.9027369297630505, False) 

Starting at 10.00, no root found, last estimate was 6.90, giving value f(6.90) = -6.366463e-12 

(-2.6575456505038764, False) 

Starting at 0.00, no root found, last estimate was -2.66, giving value f(-2.66) = 8.839758e+03 

(-8.106973924480215, False) 

Starting at -10.00, no root found, last estimate was -8.11, giving value f(-8.11) = -1.364242e-11 

(-8.106973924480215, False) 

Starting at -100.00, no root found, last estimate was -8.11, giving value f(-8.11) = -1.364242e-11 

Это выход мне нужно:

Polynomial: 5040x^0 + 1602x^1 + 1127x^2 - 214x^3 - 75x^4 + 4x^5 + 1x^6 

Derivative: 1602x^0 + 2254x^1 - 642x^2 - 300x^3 + 20x^4 + 6x^5 

(6.9027369297630505, False) 

Starting at 100.00, no root found, last estimate was 6.90,giving value f(6.90) = -2.91038e-11 

(6.9027369297630505, False) 

Starting at 10.00, no root found, last estimate was 6.90,giving value f(6.90) = -2.91038e-11 

(-2.657545650503874, False) 

Starting at 0.00, no root found, last estimate was -2.66,giving value f(-2.66) = 8.83976e+03 

(-8.106973924480215, True) 

Starting at -10.00, root found at x = -8.11, giving value f(-8.11)= 0.000000e+00 

(-8.106973924480215, True) 

Starting at -100.00, root found at x = -8.11, giving value f(-8.11)= 0.000000e+00 

Примечания: Пожалуйста, игнорируйте кортежи ошибочного вывода; Это результат моего метода Ньютона, где первым результатом является корень, а второй результат указывает, является ли он корнем или нет.

+0

Если ваши номера настолько далеки друг от друга, что вы просто работаете в пределе двойной точности. – Evert

+1

Ваша проблема на самом деле не в том фрагменте кода, который вы показали. Это просто простая оценка, но не поиск корня (например, сравнение с 0). Вероятно, вам следует сравнить с нулем * в пределах вашей доступной точности *. – Evert

ответ

1

Попробуйте это:

def evaluatePoly(poly, x_): 
    '''Evaluate the polynomial poly at x = x_ and return the result as a 
    floating-point number using Horner's rule''' 
    total= 0 
    degree =0 
    for coef in poly: 
     total += (x_**degree) * coef 
     degree += 1 
+0

Смотрите, это код, который я первоначально имел для моей программы, которая работала правильно. Это обычно оценивает полином. Тем не менее, мне нужно оценить его с помощью правила Хорнера, в котором я сталкиваюсь с моей проблемой. – schCivil

0

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

In [6]: x0=6.9027369297630505 
In [7]: evaluatePoly(poly,x0)          
Out[7]: -6.366462912410498e-12 

In [8]: evaluatePoly([abs(c) for c in poly],abs(x0)) 
Out[8]: 481315.82997756737 

Это значение представляет собой первую оценку для коэффициента увеличения ошибок с плавающей точкой, умножаются на машине эпсилон 2.22e-16 это дает ограничение на ошибки накопленная с плавающей точкой любого метода оценки из 1.07e-10, и даже оба метода оценки дают значения комфортно ниже этой границы, что указывает на то, что корень найден в возможностях формата с плавающей запятой.

Глядя на график оценки вокруг x0 один видит, что основное предположение гладкой кривой не удается при этом увеличении, а ось х проходит в прыжке, так что не лучше значение x0 не может быть найдено:

graph of polynomial evaluation close to the root