回归数续集---自阶乘

wayne

145 =1！+4！+5！
那还有没有其他的数，其十进制的各位数字的阶乘之和等于其自身的值。

mathe

容易估算，对于k位数，其各位数字阶乘之和不超过k*9!,
所以我们要求$k*9!>=10^{k-1}$,由此得出最多8位数，让后穷举吧

wayne

2# mathe
exactly right !

wayne

2# mathe

这题对于mathe来说，太水了，mathe还是解决那个 立方分数和问题吧

chyanog

这道题似曾相识，记得年初好像在哪本VB二级的书上见过，而且ProjectEuler的34题求的是所有这样的数之和。我的穷举方法，肯定算不上高明，当做是C语言的练习吧,希望看到大牛的代码出现。

1. 
   
2. #include <stdio.h>
   
3. #include <time.h>
   
4. #define max 10000000 //10^7
   
5. int main()
   
6. {        double t0=clock();
   
7.         int n,i,t,sum,fac,count=0;
   
8.         for(n=1;n<max;n++)
   
9.         {        for(sum=0,i=n;i;i/=10)
   
10.                 {        for(fac=1,t=1;t<=i%10;t++)
    
11.                         {
    
12.                                 fac*=t;
    
13.                         }
    
14.                         sum+=fac;
    
15.                 }
    
16.                 if(sum==n) printf("%d:%4d\n",++count,n);
    
17.         }
    
18.         printf("%timing: %f s\n",(clock()-t0)/CLOCKS_PER_SEC);
    
19.         return 0;
    
20. }
    

复制代码

总共找到4个：
1
2
145
40585

timing: 1.750000 s

chyanog

用Mathematica就相当简单了，只是效率差强人意：

1. 
   
2. Select[Range[10^5], Tr[IntegerDigits[#]!] == # &] // Timing
   

复制代码

{0.953, {1, 2, 145, 40585}}

chyanog

我下面的代码，参考了别人的思路，效率有了较大提升：

1. 
   
2. #include <stdio.h>
   
3. #include <math.h>
   
4. #include <time.h>
   
5. int fac(int n)
   
6. {
   
7.         return n ? n*fac(n-1) : 1;
   
8. }
   
9. int main(int argc, char *argv[])
   
10. {
    
11.         double start=clock();
    
12.         int a[10];
    
13.         for(int i=0;i<=9;i++) a[i]=fac(i);
    
14.        
    
15.         int max=10;
    
16.         for(int n=1;fac(9)*n>(int)pow(10.0,n);n++)
    
17.         max*=10;
    
18.         printf("max=%d\n",max);
    
19.        
    
20.         int n,i,sum,count=1;
    
21.         for(n=1;n<=max;n++)
    
22.         {        for(sum=0,i=n;i;i/=10)
    
23.                 sum+=a[i%10];
    
24.                 if(sum==n) printf("%d: %d\n",count++,n);
    
25.         }
    
26.         printf("timing: %f s\n",(clock()-start)/CLOCKS_PER_SEC);
    
27.         return 0;
    
28. }
    

复制代码

max=10000000
1: 1
2: 2
3: 145
4: 40585
timing: 0.476000 s

Frankenstein

用Mathematica就相当简单了，只是效率差强人意：

1. 
   
2. Select[Range[10^5], Tr[IntegerDigits[#]!] == # &] // Timing
   

复制代码

{0.953, {1, 2, 145, 40585}}
chyanog 发表于 2010-9-6 17:53

这个太吓人了，呵呵。

mathe

这个可以改进一下，改成寻找满足条件的16进制的数，也就是一个16进制数是它的16进制表示各位的阶乘之和。
如果再改成32进制数呢？

chyanog

http://www.research.att.com/~njas/sequences/A014080