一道和角谷静夫猜想相关的问题

shingyuan

最近在看EulerProject上的问题:

有一道问题是这样的:

角谷静夫猜想是说：对于一个任意的自然数，总可以通过以下步骤变为1：

$f(n) = \frac{n}{2}.$  如果 n 是偶数
$f(n) = 3 \times n + 1.$ 如果n 是奇数

对于不同的数，需要的步骤肯定不同，问题是1000000以内的数，哪个数需要的步骤最多？

我的想法比较粗糙：就是对每个数去统计.然后选取一个最长的。

对于每个数，到底需要多少次才可以到达1应该是未知的。假设统计这么多所需要的平均次数是T,那么就做了 1000000*T次循环，然后还有(2*T+1)*1000000次判断，我的代码如下;

1. int main()
   
2. {
   
3.     //open file to output
   
4.     ofstream outfile("G:/math/ans.txt");
   
5.     long ans = 0;
   
6.     long index;
   
7.     long a,b;
   
8.     for(long i = 2; i < 1000000; )
   
9.     {
   
10.              a = i;
    
11.              index = 0;
    
12.              while(true)
    
13.              {
    
14.                      if(a == 1)
    
15.                      {
    
16.                           break;
    
17.                      }
    
18.                      if(a % 2 == 0)
    
19.                      {
    
20.                           a = a/2;
    
21.                      }
    
22.                      else
    
23.                      {
    
24.                          a = a*3 + 1;
    
25.                      }
    
26.                      index++;
    
27.              }
    
28.              if(index > ans)
    
29.              {
    
30.                  ans = index;
    
31.                  b = i;
    
32.              }
    
33.             i+=2;         
    
34.     }
    
35.     if(ans < 500000)
    
36.   {
    
37.        ans = ans * 2;
    
38.       b++;
    
39. }
    
40.     outfile << "这个数是:" << ans << endl;
    
41.     outfile << "需要的次数:" << b << endl;
    
42. }

复制代码

不过貌似我的电脑根本跑不出来这个结果。是不是我的计算方法太愚蠢了。（我考虑偶数可以不用判断了。如果判断出来最大的数是<500000的奇数，则考虑这个数*2的那个偶数，就是需要步骤最多的。）我跑了很久结果也出不来

liangbch

我发现楼主的代码有一个bug。楼主显然没有考虑到整数溢出问题，在执行下面的语句时，可能会发生整数溢出而导致结果错误。例如，当n=159487，在计算到第60步时，a*3+1的值就会超过32bit,这时，如果变量类型为unsigned long, 在32bit编译器下编译，就会出现错误。我注意到楼主的代码中，数据类型为long,那么溢出的几率会更高

  a = a*3 + 1;

liangbch

　 我没有仔细调试楼主的程序。我自己写了一个程序，在PM1.7上运行，计算
100万/1000万/1亿之内的数需要的步时， 需要375/6250/93265毫秒。

　几点说明：
1.在计算过程中，我的程序使用unsigned __int64(VC编译器) 型变量，避免溢出。
2.我的程序作了一点儿优化，事先计算出1-65536之间的数所需的步骤，在计算大于65536的数所需的步骤时，一旦发现计算过程中，a的值<=65536,那么剩余的步骤可直接查表。

liangbch

本楼给出所有源代码

1. #include <stdio.h>
   
2. #include <stdlib.h>
   
3. #include <windows.h>
   
4. 
   
5. typedef unsigned __int64 UINT64;
   
6. typedef unsigned long DWORD;
   
7. typedef unsigned short WORD;
   
8. 
   
9. WORD g_steps_tab[65536];
   
10. 
    
11. int steps(DWORD n)
    
12. {
    
13.         UINT64 x=n;
    
14.         int r=0;
    
15.         while ( x >1)
    
16.         {
    
17.                 if ( x & 1)  //x 是奇数
    
18.                 {
    
19.                         if ( x>=6148914691236517205I64)
    
20.                         {
    
21.                                 printf("ERROR,n=%u,integer overflow, can not continue\n",n);
    
22.                                 return 0;
    
23.                         }
    
24.                         x=x*2+x+1;
    
25.                         r++;
    
26.                 }
    
27.                 else
    
28.                 {
    
29.                         x=x /2;
    
30.                         r++;
    
31.                 }
    
32.         }
    
33.         return r;
    
34. }
    
35. 
    
36. 
    
37. //在计算过程中,当n的值小于65536，直接查表
    
38. //tab[i] 表示 i+1 转化为1需要的步数
    
39. int steps_pro(DWORD n,WORD tab[])
    
40. {
    
41.         UINT64 x=n;
    
42.         int r=0;
    
43.        
    
44.         while ( x >1)
    
45.         {
    
46.                 if ( x <=65536)
    
47.                 {
    
48.                         return r+tab[x-1];
    
49.                 }
    
50.                 else if ( x & 1)  //x 是奇数
    
51.                 {
    
52.                         if ( x>=6148914691236517205I64)
    
53.                         {
    
54.                                 printf("ERROR,n=%u,integer overflow, can not continue\n",n);
    
55.                                 return 0;
    
56.                         }
    
57.                         x=x*2+x+1;
    
58.                         r++;
    
59.                 }
    
60.                 else
    
61.                 {
    
62.                         x=x /2;
    
63.                         r++;
    
64.                 }
    
65.         }
    
66.         return r;
    
67. }
    
68. 
    
69. void test2(DWORD to)
    
70. {
    
71.         DWORD i,max_n;
    
72.         int r,max_r;
    
73.         max_r=0;
    
74. 
    
75.         DWORD t;
    
76.         t=GetTickCount();
    
77.         for (i=1;i<=65536;i++)
    
78.         {
    
79.                 r=steps(i);
    
80.                 g_steps_tab[i-1]=r;
    
81.         }
    
82.        
    
83.         max_r=0;
    
84.         for (i=1;i<=to;i++)
    
85.         {
    
86.                 r=steps_pro(i,g_steps_tab);
    
87.                 if (r>max_r)
    
88.                 {
    
89.                         max_r=r;
    
90.                         max_n=i;
    
91.                 }
    
92.         }
    
93.         t=GetTickCount()-t;
    
94.         printf("Calc up to %u, take %d ms\n",to,t);
    
95.         printf("steps[%d]=%d\n",max_n,max_r);
    
96. }
    
97. 
    
98. 
    
99. int main(int argc, char argv[])
    
100. {
     
101.         test2(1000000);
     
102.         test2(10000000);
     
103.         test2(100000000);
     
104.         return 0;
     
105. }

复制代码

shingyuan

非常感谢~~，看了收获很大。
有个小问题：
把 x = x*3 + 1;
写成：
x=x*2+x+1;
又没有什么特别的原因？
我运行时间好像差不太多。
我是AMD 2800+ ， 512MB内存
计算1000000以内花费406MS-500MS。
判断奇数的方法+查表比我的节约非常多时间。

liangbch

是基于优化的考虑。
位运算一定比除法快。因此，
   1. 判断是否是奇数写成  "(n & 1)==1" 比 " n % 2 ==1" 快
   2.“ x=x*2+x+1;” 对通常的编译器，编译成的指令为x=(x<<1)+x+1, 不使用乘法指令。而"x=x*3+1;"对先进的编译器(如intel编译器)会编译成x=(x<<1)+x+1，但对于老的编译器，可能依然使用乘法指令来实现。乘法指令要比移位指令慢得多。

shingyuan

原来这样，好精彩，非常感谢