变态的算法题

koolcoy

昨晚想了一下这个问题，有这么一个优化方法：
将1^21 * 1, 1^21 * 2, 1^21 * 3 ... 1^21 * 21, 2^21 * 1, 2^21 * 2 ... 2^21 * 21, .... 9^21 * 21一共189个数保存起来。于是原问题抓换成从这189个数中选最多9个数出来，使得他们的和满足一定的条件。整个问题的解空间小于21^9 ，这个数字大概是5000亿，采用深搜，再根据位数为21这个条件剪一下枝，最终的解空间应该不会超过100亿。

这种算法的复杂度是O(N^9)，而不是O(9^N)，N为位数。

KBTiller

回复 81# koolcoy


    9^21不需要那么多。位数为21应该可以减去两支

koolcoy

用python实现了一下这个算法
测试结果：
[witch@localhost ~]$ for i in 21 23 24; do
> time ./flower.py $i
> done
449177399146038697307
128468643043731391252

real    1m33.091s
user    1m31.283s
sys     0m0.123s
27879694893054074471405
27907865009977052567814
35452590104031691935943
21887696841122916288858
28361281321319229463398

real    2m49.607s
user    2m47.742s
sys     0m0.164s
174088005938065293023722
188451485447897896036875
239313664430041569350093

real    3m40.527s
user    3m38.041s
sys     0m0.160s

1. 
   
2. #!/usr/bin/env python
   
3. 
   
4. import sys
   
5. 
   
6. LENGTH = int(sys.argv[1])
   
7. TABLE = tuple([tuple([x * (i**LENGTH) for x in range(1, LENGTH + 1)])
   
8.                                                 for i in range(1, 10)])
   
9. MAX = 10 ** LENGTH - 1
   
10. MIN = 10 ** (LENGTH - 1)
    
11. 
    
12. def isValidNumber(number, digit):
    
13.         tmp = [0] * 10
    
14.         for ch in str(number):
    
15.                 tmp[ord(ch) - ord('0')] += 1
    
16.         tmp.pop(0)
    
17.         return tmp == digit
    
18. 
    
19. def search(left, index, partialSum, digit):
    
20.         if index < 8 and left > 0:
    
21.                 search(left, index + 1, partialSum, digit)
    
22. 
    
23.         for count in range(1, left + 1):
    
24.                 digit[index] = count
    
25.                 value = TABLE[index][count - 1]
    
26.                 sum = partialSum + value
    
27. 
    
28.                 if sum > MAX:
    
29.                         break
    
30. 
    
31.                 if sum > MIN:
    
32.                         if isValidNumber(sum, digit):
    
33.                                 print sum
    
34.                 if index < 8 and left > count:
    
35.                         search(left - count, index + 1, sum, digit)
    
36.        
    
37.         digit[index] = 0
    
38. 
    
39. def main():
    
40.         digit = [0] * 9
    
41.         search(LENGTH, 0, 0, digit)
    
42. 
    
43. if __name__ == "__main__":
    
44.         main()
    

复制代码

koolcoy

[witch@localhost ~]$ time ./flower.py 29
23866716435523975980390369295
14607640612971980372614873089
19008174136254279995012734740
19008174136254279995012734741

real    15m38.757s
user    15m0.449s
sys     0m0.333s

believetruelove

先mark一下，今天加班太累了。

koolcoy

又用C++实现了一下这个算法，结果大大地出乎我的意料{:3_189:}

witch:~/Documents/workspace/Flower$ time ./a.out 21
449177399146038697307
128468643043731391252

real        0m9.627s
user        0m9.582s
sys        0m0.008s
witch:~/Documents/workspace/Flower$ time ./a.out 24
174088005938065293023722
188451485447897896036875
239313664430041569350093

real        0m22.139s
user        0m22.119s
sys        0m0.006s

1. 
   
2. #include <cstdio>
   
3. #include <cstring>
   
4. #include <cstdlib>
   
5. 
   
6. template<int LENGTH>
   
7. class Integer {
   
8.         typedef unsigned int U32;
   
9.         typedef unsigned long long U64;
   
10. public:
    
11.         Integer() {
    
12.                 for (int i = 0; i < LENGTH; ++i) {
    
13.                         data[i] = 0;
    
14.                 }
    
15.         }
    
16. 
    
17.         Integer(U32 that) {
    
18.                 copy(that);
    
19.         }
    
20. 
    
21.         const Integer& operator = (U32 that) {
    
22.                 copy(that);
    
23.                 return *this;
    
24.         }
    
25. 
    
26.         const Integer& operator += (const Integer& that) {
    
27.                 int carry = 0;
    
28.                 for (int i = 0; i < LENGTH; ++i) {
    
29.                         data[i] += that.data[i] + carry;
    
30.                         if (data[i] > BASE) {
    
31.                                 carry = data[i] / BASE;
    
32.                                 data[i] = data[i] % BASE;
    
33.                         } else {
    
34.                                 carry = 0;
    
35.                         }
    
36.                 }
    
37.                 return *this;
    
38.         }
    
39. 
    
40.         const Integer& operator *= (U32 that) {
    
41.                 U64 carry = 0;
    
42.                 for (int i = 0; i < LENGTH; ++i) {
    
43.                         U64 tmp = (U64)data[i] * (U64)that + carry;
    
44.                         if (tmp > BASE) {
    
45.                                 carry = tmp / BASE;
    
46.                                 data[i] = tmp % BASE;
    
47.                         } else {
    
48.                                 carry = 0;
    
49.                                 data[i] = tmp;
    
50.                         }
    
51.                 }
    
52.                 return *this;
    
53.         }
    
54. 
    
55.         Integer operator + (const Integer& that) const {
    
56.                 Integer tmp = *this;
    
57.                 return tmp += that;
    
58.         }
    
59. 
    
60.         Integer operator * (U32 that) const {
    
61.                 Integer tmp = *this;
    
62.                 return tmp *= that;
    
63.         }
    
64. 
    
65.         bool operator < (const Integer& that) const {
    
66.                 for (int i = LENGTH - 1; i >= 0; --i) {
    
67.                         if (data[i] < that.data[i]) {
    
68.                                 return true;
    
69.                         } else if (data[i] > that.data[i]) {
    
70.                                 return false;
    
71.                         }
    
72.                 }
    
73.                 return false;
    
74.         }
    
75. 
    
76.         bool operator > (const Integer& that) const {
    
77.                 return std::memcmp(data, that.data, sizeof(data)) != 0
    
78.                                 && !(*this < that);
    
79.         }
    
80. 
    
81.         void validate(char digit[9]) const {
    
82.                 char decimal[LENGTH * WIDTH + 1];
    
83.                 toString(decimal);
    
84. 
    
85.                 char actual[10] = {0};
    
86.                 for (int i = 0; i < WIDTH * LENGTH; ++i) {
    
87.                         ++actual[decimal[i] - '0'];
    
88.                 }
    
89.                 if (std::memcmp(actual + 1, digit, 9) == 0) {
    
90.                         output(decimal);
    
91.                 }
    
92.         }
    
93. 
    
94.         void toString(char* decimal) const {
    
95.                 for (int i = LENGTH - 1; i >= 0; --i) {
    
96.                         std::sprintf(decimal + WIDTH * (LENGTH - 1 - i), "%09d", data[i]);
    
97.                 }
    
98.                 decimal[LENGTH * WIDTH] = 0;
    
99.         }
    
100. 
     
101.         void show() const {
     
102.                 char decimal[LENGTH * WIDTH + 1];
     
103.                 toString(decimal);
     
104.                 printf("%s\n", decimal);
     
105.         }
     
106. 
     
107. private:
     
108.         enum {
     
109.                 BASE = 1000000000, /* This is 10^9 */
     
110.                 WIDTH = 9
     
111.         };
     
112.         U32 data[LENGTH];
     
113. 
     
114.         void copy(U32 that) {
     
115.                 data[0] = that % BASE;
     
116.                 data[1] = that / BASE;
     
117.                 for (int i = 2; i < LENGTH; ++i) {
     
118.                         data[i] = 0;
     
119.                 }
     
120.         }
     
121. 
     
122.         static void output(const char* decimal) {
     
123.                 int i = 0;
     
124.                 while (decimal[i] == '0') {
     
125.                         ++i;
     
126.                 }
     
127.                 printf("%s\n", decimal + i);
     
128.         }
     
129. };
     
130. 
     
131. typedef Integer<3> Bigint;
     
132. 
     
133. Bigint MIN;
     
134. Bigint MAX;
     
135. Bigint** TABLE;
     
136. 
     
137. void createSumTable(int length) {
     
138.         TABLE = new Bigint*[9];
     
139.         for (int i = 0; i < 9; ++i) {
     
140.                 TABLE[i] = new Bigint[length];
     
141. 
     
142.                 TABLE[i][0] = i + 1;
     
143.                 for (int j = 1; j < length; ++j) {
     
144.                         TABLE[i][0] *= (i + 1);
     
145.                 }
     
146. 
     
147.                 for (int j = 1; j < length; ++j) {
     
148.                         TABLE[i][j] = TABLE[i][0] * (j + 1);
     
149.                 }
     
150.         }
     
151. }
     
152. 
     
153. void initialize(int length) {
     
154.         MIN = 10;
     
155.         for (int i = 0; i < length - 2; ++i) {
     
156.                 MIN *= 10;
     
157.         }
     
158.         MAX = MIN * 10;
     
159. 
     
160.         createSumTable(length);
     
161. }
     
162. 
     
163. void search(int left, int index, Bigint partialSum, char* digit) {
     
164.         if (index < 8 && left > 0) {
     
165.                 search(left, index + 1, partialSum, digit);
     
166.         }
     
167.         for (int count = 1; count < left + 1; ++count) {
     
168.                 digit[index] = count;
     
169.                 Bigint sum = partialSum + TABLE[index][count - 1];
     
170. 
     
171.                 if (sum > MAX) {
     
172.                         break;
     
173.                 }
     
174.                 if (sum > MIN) {
     
175.                         sum.validate(digit);
     
176.                 }
     
177.                 if (index < 8 && left > count) {
     
178.                         search(left - count, index + 1, sum, digit);
     
179.                 }
     
180.         }
     
181.         digit[index] = 0;
     
182. }
     
183. 
     
184. int main(int argc, char** argv) {
     
185.         int length = 0;
     
186.         if (argc != 2 || sscanf(argv[1], "%d", &length) != 1 || length <= 0) {
     
187.                 std::fprintf(stderr, "Usage: %s <length>\n", argv[0]);
     
188.         }
     
189. 
     
190.         initialize(length);
     
191.         char* digit = (char*)std::calloc(length, 1);
     
192.         search(length, 0, Bigint(0), digit);
     
193.         return 0;
     
194. }
     

复制代码

koolcoy

改一行代码：typedef Integer<4> Bigint
{:3_189:}
witch:~/Documents/workspace/Flower$ time ./a.out 31
1927890457142960697580636236639
2309092682616190307509695338915
1145037275765491025924292050346

real        2m55.434s
user        2m52.737s
sys        0m0.200s

koolcoy

把sprintf调用去掉后竟然{:3_182:}{:3_182:}{:3_190:}{:3_190:}

witch:~/Documents/workspace/Flower$ time ./a.out 21
449177399146038697307
128468643043731391252

real        0m2.424s
user        0m2.030s
sys        0m0.003s
witch:~/Documents/workspace/Flower$ time ./a.out 31
1927890457142960697580636236639
2309092682616190307509695338915
1145037275765491025924292050346

real        0m26.691s
user        0m26.639s
sys        0m0.020s

1. 
   
2. #include <cstdio>
   
3. #include <cstring>
   
4. #include <cstdlib>
   
5. 
   
6. template<int LENGTH>
   
7. class Integer {
   
8.         typedef unsigned int U32;
   
9.         typedef unsigned long long U64;
   
10. public:
    
11.         Integer() {
    
12.                 for (int i = 0; i < LENGTH; ++i) {
    
13.                         data[i] = 0;
    
14.                 }
    
15.         }
    
16. 
    
17.         Integer(U32 that) {
    
18.                 copy(that);
    
19.         }
    
20. 
    
21.         const Integer& operator = (U32 that) {
    
22.                 copy(that);
    
23.                 return *this;
    
24.         }
    
25. 
    
26.         const Integer& operator += (const Integer& that) {
    
27.                 int carry = 0;
    
28.                 for (int i = 0; i < LENGTH; ++i) {
    
29.                         data[i] += that.data[i] + carry;
    
30.                         if (data[i] > BASE) {
    
31.                                 carry = data[i] / BASE;
    
32.                                 data[i] = data[i] % BASE;
    
33.                         } else {
    
34.                                 carry = 0;
    
35.                         }
    
36.                 }
    
37.                 return *this;
    
38.         }
    
39. 
    
40.         const Integer& operator *= (U32 that) {
    
41.                 U64 carry = 0;
    
42.                 for (int i = 0; i < LENGTH; ++i) {
    
43.                         U64 tmp = (U64)data[i] * (U64)that + carry;
    
44.                         if (tmp > BASE) {
    
45.                                 carry = tmp / BASE;
    
46.                                 data[i] = tmp % BASE;
    
47.                         } else {
    
48.                                 carry = 0;
    
49.                                 data[i] = tmp;
    
50.                         }
    
51.                 }
    
52.                 return *this;
    
53.         }
    
54. 
    
55.         Integer operator + (const Integer& that) const {
    
56.                 Integer tmp = *this;
    
57.                 return tmp += that;
    
58.         }
    
59. 
    
60.         Integer operator * (U32 that) const {
    
61.                 Integer tmp = *this;
    
62.                 return tmp *= that;
    
63.         }
    
64. 
    
65.         bool operator > (const Integer& that) const {
    
66.                 for (int i = LENGTH - 1; i >= 0; --i) {
    
67.                         if (data[i] > that.data[i]) {
    
68.                                 return true;
    
69.                         } else if (data[i] < that.data[i]) {
    
70.                                 return false;
    
71.                         }
    
72.                 }
    
73.                 return false;
    
74.         }
    
75. 
    
76.         bool operator < (const Integer& that) const {
    
77.                 return std::memcmp(data, that.data, sizeof(data)) != 0
    
78.                                 && !(*this > that);
    
79.         }
    
80. 
    
81.         void validate(char digit[9]) const {
    
82.                 char decimal[LENGTH * WIDTH + 1];
    
83.                 toString(decimal);
    
84. 
    
85.                 char actual[10] = {0};
    
86.                 for (int i = 0; i < WIDTH * LENGTH; ++i) {
    
87.                         ++actual[decimal[i] - '0'];
    
88.                 }
    
89.                 if (std::memcmp(actual + 1, digit, 9) == 0) {
    
90.                         output(decimal);
    
91.                 }
    
92.         }
    
93. 
    
94.         void toString(char* decimal) const {
    
95.                 for (int i = 0; i < LENGTH; ++i) {
    
96.                         convertChunk(data[LENGTH - i - 1], decimal + WIDTH * i);
    
97.                 }
    
98.                 decimal[LENGTH * WIDTH] = 0;
    
99.         }
    
100. 
     
101. private:
     
102.         enum {
     
103.                 BASE = 1000000000, /* This is 10^9 */
     
104.                 WIDTH = 9
     
105.         };
     
106.         U32 data[LENGTH];
     
107. 
     
108.         static void output(const char* decimal) {
     
109.                 int i = 0;
     
110.                 while (decimal[i] == '0') {
     
111.                         ++i;
     
112.                 }
     
113.                 printf("%s\n", decimal + i);
     
114.         }
     
115. 
     
116.         static void convertChunk(U32 digit, char* decimal) {
     
117.                 for (int i = 0; i < WIDTH; ++i) {
     
118.                         decimal[WIDTH - 1 - i] = digit % 10 + '0';
     
119.                         digit /= 10;
     
120.                 }
     
121.         }
     
122. 
     
123.         void copy(U32 that) {
     
124.                 data[0] = that % BASE;
     
125.                 data[1] = that / BASE;
     
126.                 for (int i = 2; i < LENGTH; ++i) {
     
127.                         data[i] = 0;
     
128.                 }
     
129.         }
     
130. };
     
131. 
     
132. typedef Integer<SIZE> Bigint;
     
133. 
     
134. Bigint MIN;
     
135. Bigint MAX;
     
136. Bigint** TABLE;
     
137. 
     
138. void createSumTable(int length) {
     
139.         TABLE = new Bigint*[9];
     
140.         for (int i = 0; i < 9; ++i) {
     
141.                 TABLE[i] = new Bigint[length];
     
142. 
     
143.                 TABLE[i][0] = i + 1;
     
144.                 for (int j = 1; j < length; ++j) {
     
145.                         TABLE[i][0] *= (i + 1);
     
146.                 }
     
147. 
     
148.                 for (int j = 1; j < length; ++j) {
     
149.                         TABLE[i][j] = TABLE[i][0] * (j + 1);
     
150.                 }
     
151.         }
     
152. }
     
153. 
     
154. void initialize(int length) {
     
155.         MIN = 10;
     
156.         for (int i = 0; i < length - 2; ++i) {
     
157.                 MIN *= 10;
     
158.         }
     
159.         MAX = MIN * 10;
     
160. 
     
161.         createSumTable(length);
     
162. }
     
163. 
     
164. void search(int left, int index, Bigint partialSum, char* digit) {
     
165.         if (index < 8 && left > 0) {
     
166.                 search(left, index + 1, partialSum, digit);
     
167.         }
     
168.         for (int count = 1; count < left + 1; ++count) {
     
169.                 digit[index] = count;
     
170.                 Bigint sum = partialSum + TABLE[index][count - 1];
     
171. 
     
172.                 if (sum > MAX) {
     
173.                         break;
     
174.                 }
     
175.                 if (sum > MIN) {
     
176.                         sum.validate(digit);
     
177.                 }
     
178.                 if (index < 8 && left > count) {
     
179.                         search(left - count, index + 1, sum, digit);
     
180.                 }
     
181.         }
     
182.         digit[index] = 0;
     
183. }
     
184. 
     
185. int main(int argc, char** argv) {
     
186.         int length = 0;
     
187.         if (argc != 2 || sscanf(argv[1], "%d", &length) != 1 || length <= 0) {
     
188.                 std::fprintf(stderr, "Usage: %s <length>\n", argv[0]);
     
189.         }
     
190. 
     
191.         initialize(length);
     
192.         char* digit = (char*)std::calloc(length, 1);
     
193.         search(length, 0, Bigint(0), digit);
     
194.         return 0;
     
195. }
     

复制代码

KBTiller

只要是提供组合数的解法
各种方案的解空间是一样大的

但是
for(n9 = 1 ; n9 <= 21 ; n9++)
    for(n8 = 0 ; n8 <= 21 - n9 ; n8++)
    ……
方案对所谓的“剪枝”可能更好写

koolcoy

回复 89# KBTiller


    显然不一样，按这种方法穷举的话，其复杂度是O(N^9)，如果按各个数位上的数字来枚举的话，复杂度O(9^N)，一个是多项式级的复杂度，另外一个是指数级的复杂度。复杂度跟解空间是等价的。