Выход наночастиц из-за расширения синуса маклаурина, консольные сбои - PullRequest
Новая
       7

Выход наночастиц из-за расширения синуса маклаурина, консольные сбои

0 голосов
/ 04 июня 2018

Вот мой код: 

#include <iostream>
#include <cmath>

using namespace std;

int factorial(int);

int main()
{
    for(int k = 0; k < 100000; k++)
    {
        static double sum = 0.0;
        double term;
        term = (double)pow(-1.0, k) * (double)pow(4.0, 2*k+1) / factorial(2*k+1);
        sum = sum + term;
        cout << sum << '\n';
    }
}





int factorial(int n)
{
    if(n == 0)
    {
        return 1;
    }
    return n*factorial(n-1);
}

Я просто пытаюсь вычислить значение sine(4), используя maclaurin форму расширения sine.Для каждого выхода консоли значение читается как «nan».Консоль выдает ошибку и выключается примерно через 10 секунд.Я не получаю никаких ошибок в IDE.

Ответы [ 2 ]

0 голосов
/ 05 июня 2018

Как уже говорилось другими, промежуточные результаты, которые вы получите для больших k, слишком велики, чтобы вписаться в двойную.От определенного k до pow, а также factorial будет возвращаться бесконечность.Это просто то, что происходит для очень больших двойников.И когда вы затем делите одну бесконечность на другую, вы получаете NaN.

Один из распространенных приемов для обработки слишком больших чисел - использование логарифмов для промежуточных результатов и только в конце примените экспоненциальную функцию один раз.Некоторое математическое знание логарифмов требуется здесь.Чтобы понять, что я здесь делаю, вам нужно знать exp(log(x)) == x, log(a^b) == b*log(a) и log(a/b) == log(a) - log(b).

. В вашем случае вы можете переписать 

pow(4, 2*k+1)

на 

exp((2*k+1)*log(4))

Тогда еще есть факториал.Функция lgamma может помочь с factorial(n) == gamma(n+1) и log(factorial(n)) == lgamma(n+1).Короче говоря, lgamma дает вам журнал факториала без огромных промежуточных результатов.

Подводя итог, замените 

pow(4, 2*k+1) / factorial(2*k+1)

на 

exp((2*k+1)*log(4) - lgamma(2*k+2))

Это должно помочь вамс вашими NaNs.Кроме того, это должно повысить производительность, поскольку lgamma работает в O(1), тогда как ваш factorial находится в O(k).

Обратите внимание, однако, что я все еще очень мало уверен в том, что ваш результат будет численно точным.Точность двойного кода по-прежнему ограничена примерно 16 десятичными цифрами.Ваши 100000 итераций, скорее всего, бесполезны, возможно, даже вредны.

0 голосов
/ 04 июня 2018

В вашем подходе много проблем.

Ваша факторная функция не может вернуть int.Возвращаемое значение будет слишком большим, очень быстрым.

Использование pow(-1, value) для получения чередующегося положительного / отрицательного значения очень неэффективно и довольно быстро даст неправильное значение.Вы должны выбрать 1,0 или -1,0 в зависимости от четности k.

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

Алгоритмически, если вы собираетесь поднять 4 до степени 2k + 1, а затем разделить на (2k + 1) !, вы должны оставить обасписок факторов (4, 4, 4, 4 ...) и (2,3,4,5,6,7,8,9, ....) и упрощают обе стороны.На числителях и знаменателях одновременно можно удалить четыре четверки.

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

...