Вы спрашивали, есть ли 10!слишком велик, но если бы вы наблюдали i в точке переполнения значения, вы увидите, что i == 12.Тогда выражение для temp будет:
3 25 / 25 x 12!
Подвыражение 25 x 12!имеет значение 11975040000, для которого потребуется 34 бита, если оно является целочисленным выражением, которое будет, если fact() является целочисленной функцией.Если это 32-разрядная целочисленная функция, она будет переполнена.
Если fact() изменить на 64-разрядный тип без знака, он получит дальнейшее развитие, но это хорошо только для 19!и подвыражение делителя переполняется в i=18 (29 x 18!) перед сходимостью.
Он будет сходиться для i = 3, используя следующую реализацию факториальной функции:
double fact(int n)
{
if (n >= 1)
return n * fact(n - 1);
else
return 1;
}
С результатом 0,886207.Для произвольно больших значений x вам понадобится математическая библиотека произвольной точности или «bignum», а не встроенные типы.
Вдохновленный комментарием Ли Крокера, я пробовал различные выражения (проверено на GDB Online )
1) Исходя из моих double fact( int ) и менее ваших посторонних скобок;
double e = (2*i)+1 ;
temp = pow(x,e) / (e * fact(i));
Результаты для x = 1 до 9:
x: 1.0 s(x): 0.746824
x: 2.0 s(x): 0.882081
x: 3.0 s(x): 0.886207
x: 4.0 s(x): 0.886227
x: 5.0 s(x): 0.886227
x: 6.0 s(x): 0.889677
x: 7.0 s(x): 2453.707036
x: 8.0 s(x): 1055556860.115990
x: 9.0 s(x): -nan
2) Использование tgamma():
temp = pow(x,e) / (e * tgamma(i+1)) ;
Результат:
x: 1.0 s(x): 0.746824
x: 2.0 s(x): 0.882081
x: 3.0 s(x): 0.886207
x: 4.0 s(x): 0.886227
x: 5.0 s(x): 0.886226
x: 6.0 s(x): 0.882858
x: 7.0 s(x): -83.077451
x: 8.0 s(x): 3729162898.780624
3) Использование lgamma():
temp = pow(x,e) / exp(log(e) + lgamma(i+1)) ;
Результат:
x: 1.0 s(x): 0.746824
x: 2.0 s(x): 0.882081
x: 3.0 s(x): 0.886207
x: 4.0 s(x): 0.886227
x: 5.0 s(x): 0.886236
x: 6.0 s(x): 1.904565
x: 7.0 s(x): -322346.422518
x: 8.0 s(x): 5416428382.734000
x: 9.0 s(x): -nan
Результаты отличаются для x = 6.0, но все версии разбиты примерно там.Что лучше, я не ясно - fact() Я подозреваю.Стандартные библиотечные функции не предлагают никаких существенных улучшений.
4) Добавлено выражение Ли для проверки (если я правильно понял его ответ):
temp = exp(e * log(x) - log(e) - lgamma(i+1)) ;
Результат:
x: 1.0 s(x): 0.746824
x: 2.0 s(x): 0.882081
x: 3.0 s(x): 0.886207
x: 4.0 s(x): 0.886227
x: 5.0 s(x): 0.886245
x: 6.0 s(x): 2.892730
x: 7.0 s(x): 96596.597282
x: 8.0 s(x): -720878097610.696167
x: 9.0 s(x): -18109761159507625984.000000
Для справки мой тестовый код:
double s( double x )
{
double sx = x ;
double temp = 1.0 ;
for( int i = 1; temp > 1e-6; i++ )
{
double e = (2*i)+1 ;
// Uncomment one of the following:
//temp = pow(x,e) / (e * fact(i));
//temp = pow(x,e) / (e * tgamma(i+1)) ;
//temp = pow(x,e) / exp(log(e) + lgamma(i+1)) ;
//temp = exp(e * log(x) - log(e) - lgamma(i+1)) ;
sx += (i % 2) == 0 ? temp : -temp ;
}
return sx ;
}
int main()
{
for( double x = 1.0; x < 10.0; x += 1.0 )
{
printf( "x: %.1f s(x): %lf\n", x, s(x) ) ;
}
return 0;
}