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

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>
5. typedef unsigned __int64 UINT64;
6. typedef unsigned long DWORD;
7. typedef unsigned short WORD;
9. WORD g_steps_tab[65536];
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. }
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;
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. }
69. void test2(DWORD to)
70. {
71. DWORD i,max_n;
72. int r,max_r;
73. max_r=0;
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. }
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. }
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

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