平方数数字和

shshsh_0510

看mathe的方法，对于很多的x，应该都可以构造01向量来达到。 那么，接下来的困难得多的问题是，对于小于N的自然数x=a的大概有多少呢？

mathe

这个问题比较难。可靠的方法很难想到。 不过我们可以假设平方数各位数都是均匀分布的，给出一个近似统计数据来（计算公式估计还是有点麻烦）。 而更加精确的模型可以对平方数的首位和末位的分布各个更加精确一些的统计模型，比如末位应该是$1/5$的概率取1,4,9,6,然后$1/10$的概率取0和5

gxqcn

刚才把 36# 算法再次优化， 完全消除了除法及64位的运算， 算出 X=130，仅需 11.243 s；遍历18位以内的平方数，仅需 126.074 s 大家看看还有哪里可以优化的？

mathe

函数SumOfDigits可以提速的,比如我们事先计算1~9999所有数字的sumofdigit,并且保存 然后对于每个整数,每四位的对应和再累加,可以减少一些乘除次数.

gxqcn

类似想法我在 38# 曾提到过。

无心人

我想 折半算法更好吧 由于这个求数位和是并行操作 所以考虑如果采取精细的汇编 假设考虑到32位十进制数字 则首先除以10^16 分开，分开的部分 分别再除以10^8 则只需要 一次除以10^16 二次除以10^8 四次除以10^4

无心人

更黑的方法在下面 如果一开始就存储 1-9999内的乘法 共5*10^7个结果 则可以考虑整体的数字 采取10000进制 呵呵

无心人

但上述方法只能对计算超过2^64以上的 数字的平方和才有好的效果

gxqcn

下面的代码可以遍历全部的64bit内的完全平方数，且效率非常高：

1. #include <stdlib.h>
2. #include <stdio.h>
4. #include "../../../HugeCalc_API/CppAPI/Include/HugeCalc.h" // 公共接口
6. #pragma message( "automatic link to ../../../HugeCalc_API/CppAPI/Lib/HugeCalc.lib" )
7. #pragma comment( lib, "../../../HugeCalc_API/CppAPI/Lib/HugeCalc.lib" )
10. #define TEN_POW2 100UL
11. #define TEN_POW4 10000UL
13. int main( void )
14. {
15. const UINT32 cycleDelta[ 4 ] = { 1, 3, 3, 2 }; // 0, 1, 4, 7
16. UINT32 table[ TEN_POW4 ];
17. UINT32 mark[ 9*4*5 + 1 ];
18. UINT32 value[ 5 ], delta[ 5 ];
19. UINT32 * p = table;
21. UINT32 i, j, k, m, s;
22. UINT32 i2, i3, i4, s2, s3;
24. HugeCalc::EnableTimer();
25. HugeCalc::ResetTimer();
27. for ( i4 = 0; i4 < 10; ++i4 )
28. {
29. for ( i3 = 0; i3 < 10; ++i3 )
30. {
31. s3 = i4 + i3;
32. for ( i2 = 0; i2 < 10; ++i2 )
33. {
34. s2 = s3 + i2;
35. for ( i = 0; i < 10; ++i )
36. {
37. *p++ = s2 + i;
38. }
39. }
40. }
41. }
43. memset( mark, 0, sizeof( mark ));
44. memset( value, 0, sizeof( value ));
45. memset( delta, 0, sizeof( delta ));
47. k = 1;
48. j = TEN_POW2;
49. value[ 0 ] = 0;
50. delta[ 0 ] = 1;
51. m = 1;
52. for ( ; ; )
53. {
54. s = 0;
55. for ( i = 0; i < m; ++i )
56. {
57. value[ i ] += delta[ i ];
58. if ( value[ i ] >= TEN_POW4 )
59. {
60. value[ i ] -= TEN_POW4;
61. ++value[ i + 1 ];
62. }
64. s += table[ value[ i ] ];
65. }
66. if ( 0 == mark[ s ] )
67. {
68. mark[ s ] = k;
69. }
72. if ( j == ++k )
73. {
74. if ( 0 == k )
75. {
76. break;
77. }
79. j *= ( 5 == ++m ) ? 0 : TEN_POW2;
80. }
82. delta[ 0 ] += 2;
84. i = 0;
85. while ( delta[ i ] >= TEN_POW4 )
86. {
87. delta[ i ] -= TEN_POW4;
88. ++delta[ ++i ];
89. }
90. }
92. HugeCalc::EnableTimer( FALSE );
95. s = 0;
96. i = 0;
97. while ( s < 9*4*5 )
98. {
99. s += cycleDelta[ i ];
100. if ( 4 == ++i )
101. {
102. i = 0;
103. }
105. if ( 0 != mark[ s ] )
106. {
107. if ( mark[ s ] <= 0xFFFF )
108. {
109. printf( "X = %u\tS = %u\n", s, mark[ s ] * mark[ s ] );
110. }
111. else
112. {
113. printf( "X = %u\tS = %I64u\n", s, UInt32x32To64( mark[ s ], mark[ s ] ));
114. }
115. }
116. else
117. {
118. printf( "X = %u\n", s );
119. }
120. }
122. printf( HugeCalc::GetTimerStr( FT_DOT06SEC_s ) );
123. printf( "\n" );
125. return 0;
126. }

复制代码

备注：代码中，仅调用了 HugeCalc 的高精度计时部分。

gxqcn

X = 1 S = 1 X = 4 S = 4 X = 7 S = 16 X = 9 S = 9 X = 10 S = 64 X = 13 S = 49 X = 16 S = 169 X = 18 S = 576 X = 19 S = 289 X = 22 S = 1849 X = 25 S = 4489 X = 27 S = 3969 X = 28 S = 17956 X = 31 S = 6889 X = 34 S = 27889 X = 36 S = 69696 X = 37 S = 98596 X = 40 S = 97969 X = 43 S = 499849 X = 45 S = 1887876 X = 46 S = 698896 X = 49 S = 2778889 X = 52 S = 4999696 X = 54 S = 9696996 X = 55 S = 19998784 X = 58 S = 46689889 X = 61 S = 66699889 X = 63 S = 79869969 X = 64 S = 277788889 X = 67 S = 478996996 X = 70 S = 876988996 X = 72 S = 3679999569 X = 73 S = 1749999889 X = 76 S = 5599977889 X = 79 S = 7998976969 X = 81 S = 8998988769 X = 82 S = 17999978896 X = 85 S = 36799899889 X = 88 S = 88998998929 X = 90 S = 297889998849 X = 91 S = 299879997769 X = 94 S = 897977978689 X = 97 S = 975979998889 X = 99 S = 3957779999889 X = 100 S = 2699997789889 X = 103 S = 9879498789889 X = 106 S = 9998768898889 X = 108 S = 29998985899689 X = 109 S = 85986989688889 X = 112 S = 97888999968769 X = 115 S = 386999898769969 X = 117 S = 429998989997889 X = 118 S = 578889999977689 X = 121 S = 898999897988929 X = 124 S = 1959999889996996 X = 126 S = 6788999798879769 X = 127 S = 3699998989898689 X = 130 S = 9895699989899689 X = 133 S = 38896878989988889 X = 135 S = 67699789959899889 X = 136 S = 38999699989995889 X = 139 S = 188997899869998769 X = 142 S = 279869897899999969 X = 144 S = 498999778899898896 X = 145 S = 1877896979979898969 X = 148 S = 989879999979599689 X = 151 S = 6979497898999879969 X = 153 S = 5899989587897999889 X = 154 S = 8899988895999696889 X = 157 X = 160 X = 162 X = 163 X = 166 X = 169 X = 171 X = 172 X = 175 X = 178 X = 180 103.824245 s Press any key to continue

相对于 37#，不仅把 X=145 这个空洞给填补了，还多了更大的三组数据。 发过帖后，才发现本算法已实现了 无心人 在 47# 中提到的 10000 进制进行优化。 注意，我的算法中完全没有64bit长整型，不含除法，连乘法都极少。