2016-01-18 3 views
0

Я пытаюсь использовать интеграцию Монте-Карло, чтобы приблизить область под заданным графом, чтобы вычислить ее площадь. Чтобы это работало, мне нужно, чтобы вычисленные y_min и y_max были точными. Итак, в качестве примера я буду использовать график sin(x) от 0 до pi. Для того, чтобы найти y_min и y_max У меня есть следующие функции:Избегайте неточности при использовании чисел с плавающей запятой?

def y_range(f, x_min, x_max, n=100): 
    # Step size 
    h = float((x_max - x_min))/n 

    # Calculate y for n points between x_min and x_max 
    y = [f(x * h) for x in range(0, n + 1)] 

    # Get minimum and maximum y 
    y_max = max(y) 
    y_min = min(y) 

    return y_min, y_max 

Печать y_min и y_max дает:

y_max = 1.0 
y_min = -3.21624529935e-16 

Я знаю y_min должна равняться 0.0, так как я могу исправить эту неточность?

+4

Это * действительно близко * к нулю ... – jonrsharpe

ответ

0

Я бы не использовал float, когда точность важна. Предлагаю вам посмотреть decimal.Decimal. Предполагая, что вы измените x_min и x_max в десятичном, вы должны удалить только что float на линии 3.

This пост может помочь объяснить, где неточность исходит от:

+0

'decimal.Decimal' не поможет вам точно представлять pi. – Sneftel

+0

@Sneftel Хорошая точка. Возможно, мой ответ здесь не применим. –

0

пи нерационально, поэтому вы должны принять что он может быть представлен только с определенной точностью, независимо от того, какой вид представления вы используете. Если вы используете float, то точность будет в порядке 10**(-16), как указано here. Самый простой путь вперед - округлить результат до 15 десятичных цифр с помощью round(value,15). Если вам нужна более высокую точность, чем 15 десятичных цифр, вы можете действительно смотреть на другие типы для представления переменных, таких как десятичного

0

Основная проблема в том, что max не может быть получена в результате min + 100*h для любого h. Существует несколько потенциальных причин, но наиболее простым является то, что количество шагов между ними не делится на 100.

Как точно сделать это лучше зависит от того, насколько вы будете осторожны. Большие две вещи, которые вам нужны, - это интерполяция между этими двумя значениями (а не основанная на начальной точке и шаге) и точное выполнение интерполяции.

Следующий код будет производить надежные результаты:

def interp_at_step(a, b, i, n): 
    # separate calculation of alpha and beta to avoid catastrophic cancellation 
    alpha = (n-i)/n 
    beta = i/n 
    return a*alpha + b*beta 

def y_range(f, x_min, x_max, n=100): 
    # Calculate y for n points between x_min and x_max 
    y = [f(interp_at_step(x_min, x_max, x, n)) for x in range(0, n + 1)] 

    # Get minimum and maximum y 
    y_max = max(y) 
    y_min = min(y) 

    return y_min, y_max 

Конечно, как уже упоминалось Denys, пи не может быть точно представлено. Таким образом, это облегчит ошибки из интерполяции, но не обязательно из самих операндов.