воскресенье, 25 сентября 2011 г.

Вычисление функции ошибок erf() в C/C++

Как оказалось, C/C++ не содержат стандартных функций для вычисления функции ошибок. Есть три выхода из данной ситуации:
  1. Воспользоваться сторонними библиотеками типа Boost, которые реализуют данную функцию.
  2. Стандарт языка C99 содержит реализацию данной функции в библиотеке cmath. Поэтому, если ваш компилятор реализует данный стандарт, то подключайте данную библиотеку и пользуйтесь функцией erf(). В частности, компилятор GCC реализует эту функцию. Для пользователей Windows могу посоветовать среду разработки Code::Blocks, которая содержит в себе данный компилятор.
  3. Если всё-таки ваш компилятор не реализует, то можно воспользоваться следующей реализацией:

double erf(double x)
{
    /* 
        erf(z) = 2/sqrt(pi) * Integral(0..x) exp( -t^2) dt
        erf(0.01) = 0.0112834772 erf(3.7) = 0.9999998325
        Abramowitz/Stegun: p299, |erf(z)-erf| <= 1.5*10^(-7)
    */

    double y = 1.0 / ( 1.0 + 0.3275911 * x);
    return 1 - (((((
        + 1.061405429  * y
        - 1.453152027) * y
        + 1.421413741) * y
        - 0.284496736) * y
        + 0.254829592) * y)
        * exp (-x * x);
}

Приведённая выше аппроксимация действительна только для x>=0, а для остальной области аргумента значения функции erf получаются аналогично из соображений нечетности:
double erf(double x)
{
    double y = 1.0 / ( 1.0 + 0.3275911 * (x >= 0 ? x : -x));
    y = 1 - (((((
        + 1.061405429  * y
        - 1.453152027) * y
        + 1.421413741) * y
        - 0.284496736) * y
        + 0.254829592) * y)
        * exp (-x * x);
    return (x >= 0 ? y : -y);
}

4 комментария:

  1. Здесь опечатка: используется аппроксимация 7.1.26 на той же странице.

    ОтветитьУдалить
  2. я бы ещё специально для таких тупиц как я подписал бы что это апроксимация только для x>=0, а для остальной области аналогично из соображений нечетности.

    ОтветитьУдалить
    Ответы
    1. :)

      Спасибо за замечание! Добавил вторую реализацию с учётом вашего замечания.

      Удалить