вычисление функции в Matlab с очень маленькими значениями

Проблема в том, что exp(x) примерно равно 1+x, но оценивается как 1 из-за того, что 1 неотличимо от 1+x в представлении с плавающей запятой. Существует функция MATLAB expm1(x) (которая exp(x)-1 реализована для малых x), которая позволяет избежать этой проблемы и хорошо работает для малых аргументов:

>> x=1e-100;expm1(x)/x
ans =
     1

Мне пришлось попробовать мой инструмент LIMEST. Как и любой адаптивный инструмент, его можно обмануть, но обычно он довольно хорош.

fun = @(x) (exp(x) - 1)./x;

Как видите, у веселья есть проблемы на нуле.

fun(0)
ans =
   NaN

хотя если мы оцениваем веселье около нуля, мы видим, что оно близко к 1.

format long g
fun(1e-5)
ans =
          1.00000500000696

LIMEST преуспевает и даже может предоставить оценку ошибки предела.

[lim,err] = limest(fun,0,'methodorder',3)

lim =
     1

err =
      2.50668568491927e-15

LIMEST использует последовательность полиномиальных аппроксимаций в сочетании с адаптивной экстраполяцией Ричардсона для получения как предельной оценки, так и меры неопределенности этого предела.

Так какую проблему вы видите? Ошибка, которую вы видели, — это простая ошибка отмены вычитания. Таким образом, посмотрите на значение

exp(1e-20)
ans =
     1

Даже с форматом long g значение двойной точности exp(1e-20) просто слишком близко к 1, поэтому, когда мы вычитаем из 1, результатом будет точный ноль. Разделите это на любое ненулевое значение, и мы получим ноль. Конечно, когда x на самом деле равен нулю, у нас есть условие 0/0, поэтому, когда я попробовал это, получилось NaN.

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

DefaultNumberOfDigits 64
exp(hpf('1e-20'))
ans =
1.000000000000000000010000000000000000000050000000000000000000166

Обратите внимание, что когда мы вычитаем 1, разница между 1 и экспоненциальным значением меньше, чем eps(1), поэтому MATLAB должен выдать нулевое значение.

exp(hpf('1e-20')) - 1
ans =
1.000000000000000000005000000000000000000016666666666670000000000e-20

Незаданный вопрос заключается в том, как бы я выбрал точную генерацию этой функции в MATLAB. Ясно, что вы не хотите использовать грубую силу и определять удовольствие, как я, поскольку вы теряете большую часть точности. Я бы, вероятно, просто расширил ряд Тейлора в ограниченной области вокруг нуля и использовал fun, как указано выше, для x, значительно отличающегося от нуля.

Я думаю, это связано с точностью ваших чисел. Короче говоря, точность по умолчанию для MATLAB составляет 5 цифр. Вы можете расширить его до 15 с помощью format long. Прочтите эту статью для получения дополнительной информации о точности MATLAB.