double MyClass::dx = ?????;
double MyClass::f(double x)
{
return 3.0*x*x*x - 2.0*x*x + x - 5.0;
}
double MyClass::fp(double x) // derivative of f(x), that is f'(x)
{
return (f(x + dx) - f(x)) / dx;
}
При использовании метода конечных разностей для деривации очень важно выбрать оптимальный dx
значение. Математически, dx
должно быть как можно меньше. Тем не менее, я не уверен, что это правильный выбор — выбрать наименьшее положительное число двойной точности (т. Е. 2,2250738585072014 x 10).-308).
Существует ли оптимальный числовой интервал или точное значение для выбора dx
чтобы ошибка вычисления была как можно меньше?
(Я использую 64-битный компилятор. Я буду запускать свою программу на процессоре Intel i5.)
Выбор наименьшего возможного значения почти наверняка неверен: если dx
это было наименьшее число, то f(x + dx)
было бы именно так равно f(x)
из-за округления.
Таким образом, у вас есть компромисс: выберите dx
слишком мал, и вы теряете точность из-за ошибок округления. Выберите его слишком большим, и ваш результат будет неточным из-за изменений в производной Икс изменения.
Чтобы судить о числовых ошибках, рассмотрим (f(x + dx) - f(x))/f(x)
1 математически. Числитель обозначает разницу, которую вы хотите вычислить, но знаменатель обозначает величину чисел, с которыми вы имеете дело. Если эта доля составляет около 2—К, тогда вы можете ожидать примерно К немного точности в вашем результате.
Если вы знаете свою функцию, вы можете вычислить, какую ошибку вы получите, выбрав dx
слишком большой. Затем вы можете связать вещи так, чтобы ошибка, возникшая при этом, была примерно такой же, как ошибка, возникшая при округлении. Но если вы знаете функцию, вам может быть лучше, если вы предоставите функцию, которая непосредственно вычисляет производную, как в вашем примере с полигональной f
,
Раздел википедии тот погорский указал предлагает значение sqrt (ε)Икс, или примерно 1.5e-8 * x
, Без каких-либо более подробных знаний о функции такое практическое правило обеспечит разумное значение по умолчанию. Также обратите внимание, что тот же раздел предлагает не делить на dx
, но вместо этого (x + dx) - x
, так как это принимает ошибки округления, возникшие в результате вычислений x + dx
в учетную запись. Но я думаю, что вся статья полна предложений, которые вы могли бы использовать.
1 Эта формула действительно должна делиться на f(x)
не dx
хотя бывший редактор думал иначе. Я пытаюсь сравнить количество значащих битов, оставшихся после деления, а не наклон касательной.
Почему бы просто не использовать Power Power Rule для получения производной, вы получите точный ответ:
f(x) = 3x^3 - 2x^2 + x - 5
f'(x) = 9x^2 - 4x + 1
Следовательно:
f(x) = 3.0 * x * x * x - 2.0 * x * x + x - 5.0
fp(x) = 9.0 * x * x - 4.0 * x + 1.0