阶乘和圆周率

mathe

乘法运算精度损失不高的.
而使用Stirling公式计算n!,同样有效精度损失也会很低的,最主要问题来自公式本身,可能需要多取几项余项

mathe

由于
Prelude> log(3.1415926)/log(10)
0.4971498652858686
Prelude> log(3.1415927)/log(10)
0.4971498791098912
所以我们只需要搜索log(n!)/log(10)的尾数在
0.497149865和0.497149880之间的所有n就可以了。
对于搜索这个范围内的所有n,如果采用Stiring公式，那么double类型数据足够了。然后我们对它们逐一进行验算就可以了。

无心人

关键是Striling公式其误差公式是？

无心人

呵呵
用Haskell解这个问题可就简单了
不过Haskell的高精度浮点库要另外安装·

无心人

http://zh.wikipedia.org/wiki/%E6%96%AF%E7%89%B9%E6%9E%97%E5%85%AC%E5%BC%8F

上面的链接有更详细的公式样式和误差分析

mathe

对于比较大的n,使用stirling级数精度可以达到非常之高:
$ln\Gamma(z)=1/2ln(2\pi)+(z-1/2)ln(z)-z+\sum_{n=1}^{infty}\frac{B_{2n}}{2n(2n-1)z^{2n-1}}$
$=1/2ln(2\pi)+(z-1/2)ln(z)-z+1/{12z}-1/{360z^3}+1/{1260z^5}-...$
其中$B_{2n}$为贝努利数。

liangbch

因为double型数可精确到15-16位有小数字，如果斯特林公式是精确的，使用斯特林公式来计算log10(n!) 可精确到15-16位有效数字，而这15-16位有效数字又分成整数部分和小数部分，现在的要求尾数（小数部分必须达到8位），那么整数部分至多就不能超过7位了，所有当n!大于$10^7$时，就无法使用double型来估算了，必须采用更高的精度。

    而 斯特林公式（一阶形式 ）的精度有这样的特点，当n很大时，精度较高，当n较小时，精度较低。而n很大时，又会挤占小数部分的精度。故处于一个两难的境地。所以我们只能采用斯特林公式的高阶形式，在n不太大的时候来求解。当n很大时，double 型数就不行了，必须采用更高的精度。

mathe

比如现在$n>10^5$,那么取到$1/{12z}$这一项以后，误差将小于$1/{360z^3}<10^-17$,其中$z=n+1$,精度已经超过double可以表达的范围了。

mathe

原帖由 liangbch 于 2008-12-5 13:05 发表
因为double型数可精确到15-16位有小数字，如果斯特林是精确的，使用斯特林公式来计算log10(n!) 可精确到15-16位有效数字，而这15-16位有效数字又分成整数部分和小数部分，现在的要求尾数（小数部分必须达到8位），那 ...

其实主要问题在于计算${z\lg(z)$的小数部分。这个我觉得可以在算出$z\lg(z)$一个近似值u以后,可以通过迭代
$d=d-\frac{10^{d/z}-z*10^{-u/z}}{ln(10)10^{d/z}-10^{-u/z}}$
来求出d的比较精确的值（取$d_0=0$),其中对于分子，应该需要进行泰勒展开

无心人

开始记录结果
d = 3, n = 9
d = 31, n = 62
d = 314, n = 62
d = 3141, n = 10044
d =31415, n = 50583
d = 314159, n = 1490717
d = 3141592, n =   5573998  开始数字是3141592381
d = 31415926, n = 65630447 开始数字是3141592602
d = 314159265, n = 688395641 开始数字3141592651957527
d = 3141592653, n = 5777940569 起始数字314159265301
d = 31415926535, n = 77773146302, 开始 数字314159265353488
d = 314159265358,  n = 1154318938997开始数字3141592653584138
d = 3141592653589, n = 1544607046599开始数字31415926535899498