help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Bug or numerical rounding/truncation errors in factorial.m function?


From: José Luis García Pallero
Subject: Bug or numerical rounding/truncation errors in factorial.m function?
Date: Sun, 3 Oct 2010 13:35:37 +0200

Hello,
I'm using Octave 3.2.4 in a Debian Sid box. I'm trying to compute some
factorials and here is the results for some values:

factorial(16)
ans =  20922789888000
prod(1:16)
ans =  20922789888000

factorial(17)
ans =  355687428096001 *****************different
prod(1:17)
ans =  355687428096000 *****************different

factorial(18)
ans =  6402373705727994 *****************different
prod(1:18)
ans =  6402373705728000 *****************different

factorial(19)
ans =  121645100408832080 *****************different
prod(1:19)
ans =  121645100408832000 *****************different

factorial(20)
ans =  2432902008176640000
prod(1:20)
ans =  2432902008176640000

17, 18 and 19 are sufficient small for obtain the correct result via
direct multiplication. In the code of factorial.m function the
computation is performed as

x = round (gamma (n+1));

Is possible that gamma function has a bug or is only a problem with
ieee arithmetic? I suggest for the computation of factorials for small
numbers a behavior similar to the used in GSL, in the file
specfun/gamma.c (sorry, I can't paste a link because savanna don't
permit surfing the src). In this file, the factorials for numbers from
0 to 170 are stored explicitly (hard coded) for speed.

-- 
*****************************************
José Luis García Pallero
address@hidden
(o<
/ / \
V_/_
Use Debian GNU/Linux and enjoy!
*****************************************



reply via email to

[Prev in Thread] Current Thread [Next in Thread]