Вот мой код:
#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.
Есть несколько проблем с вашим подходом.
Ваша факторная функция не может вернуть 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, которые вы установили, без специального кода.
Как уже говорили другие, промежуточные результаты вы получите за большие 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 итераций, скорее всего, бесполезны, возможно, даже вредны.