阶乘有多少位？
一个数的阶乘有多少位，有好的算法吗比如四百万的阶乘有 24671066位 本帖最后由 KeyTo9_Fans 于 2010-1-20 00:11 编辑
这是什么公式来着？
$n! = \sqrt(2\pi(n+1/6+1/(72n))) * (n/e)^n$
其中$1/6$和$1/(72n)$这两项是我自己凑的。
据测试，多了那两项可以修正结果，以至于有效数字比原来多了2倍。
两边取对数，得
$log(n!) = log(2\pi(n+1/6+1/(72n)))/2 + n*(log(n)-log(e))$
把$n=4000000$代进去算，得到
$log(4000000!) = 24671065.737818781141945821796116$
大概这几位数字
$24671065.7378187811419458217961$
是有效的吧，后面的$...16$可能就不精确了。
即便如此，我们已经可以得到精确值的位数，以及精确值的前面若干位有效数字。 Stirling's formula 斯特林级数公式:
http://mathworld.wolfram.com/StirlingsSeries.html 哦，谢谢各位了 这是什么公式来着？
$n! = \sqrt(2\pi(n+1/6+1/(72n))) * (n/e)^n$
其中$1/6$和$1/(72n)$这两项是我自己凑的。
$log ...
KeyTo9_Fans 发表于 2010-1-20 00:01 http://bbs.emath.ac.cn/images/common/back.gif
这个KeyTo9_Fans 公式虽然在逼近n!上还可以有改进的空间，
但用来求阶乘的位数，恰到好处，已经是简洁到不能再简洁的了！！！ 本帖最后由 wayne 于 2010-1-20 13:35 编辑
哦，用来求大于1的数的阶乘的位数，
Stirling's formula 已经足够了
$\lceil n(\log _{10}n-\log _{10}\e)+\frac{1}{2}\log _{10}(2\pi n)\rceil $
而逼近n!,更好的是：
$n! = \sqrt(2\pi(n+1/6+1/(36.4+70.3n))) * (n/e)^n$ 本帖最后由 KeyTo9_Fans 于 2010-1-20 19:48 编辑
关于根号里第3项的分母，楼上是怎么得到$36.4+70.3n$的？
我通过计算1000000000的阶乘的准确值(仅保留100位有效数字)，然后倒过来推根号里面的多项式，可以得到
$n! =sqrt(2\pi p)*(n/e)^n$
其中
$p=n+1/6+a/n+b/n^2+c/n^3+d/n^4+e/n^5+f/n^6+g/n^7$
$a = 1/72$
$b = -31/6480$
$c = -139/155520$
$d = 9871/6531840$
$e = 324179/1175731200$
$f = -8225671/7054387200$
$g = -69685339/338610585600$
到这里，已经把100位有效数字全部用完了，无法再继续考察下一项了。
要得到更多的项，需要计算比1000000000更大的阶乘的准确值，并保留更多的有效数字。
把$n=16$代进去算，得到的结果为
$20922789887999.748608502448646388$
四舍五入到整数，与真实值
$16! =20922789888000$
完全相同。
因为根号里的多项式有9项，所以上述公式的大致趋势是$n$每增大10倍，$n!$的有效数字增加9个。
本帖最后由 wayne 于 2010-1-20 20:54 编辑
用C打印大数的阶乘，终于完成了，刚出锅，还是热的，高手们看看还有没有可以改进的地方
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
int main(){
int i,j,*f,tmp,c=0;
long int n,bits;
const double PI=2*asin(1.0),E=exp(1.0);
scanf("%ld",&n);
bits=(long)ceil(n*(log10(n)-log10(E))+log10(2*PI*n)/2);
printf("there are %ld digits\n%ld!=",bits,n);
f=(int*)calloc(bits,sizeof(int));
f=1;
for(i=2;i<=n;i++){
for(j=0;j<bits;j++){
tmp=f*i+c;
c=tmp/10;
f=tmp%10;
} }
for(i=bits-1;i>=0;i--) printf("%d",f);
printf("\n");
return 0;
}
页:
[1]
2