1 /* enough.c -- determine the maximum size of inflate's Huffman code tables over
2 * all possible valid and complete Huffman codes, subject to a length limit.
3 * Copyright (C) 2007, 2008 Mark Adler
4 * Version 1.3 17 February 2008 Mark Adler
5 */
6
7 /* Version history:
8 1.0 3 Jan 2007 First version (derived from codecount.c version 1.4)
9 1.1 4 Jan 2007 Use faster incremental table usage computation
10 Prune examine() search on previously visited states
11 1.2 5 Jan 2007 Comments clean up
12 As inflate does, decrease root for short codes
13 Refuse cases where inflate would increase root
14 1.3 17 Feb 2008 Add argument for initial root table size
15 Fix bug for initial root table size == max - 1
16 Use a macro to compute the history index
17 */
18
19 /*
20 Examine all possible Huffman codes for a given number of symbols and a
21 maximum code length in bits to determine the maximum table size for zilb's
22 inflate. Only complete Huffman codes are counted.
23
24 Two codes are considered distinct if the vectors of the number of codes per
25 length are not identical. So permutations of the symbol assignments result
26 in the same code for the counting, as do permutations of the assignments of
27 the bit values to the codes (i.e. only canonical codes are counted).
28
29 We build a code from shorter to longer lengths, determining how many symbols
30 are coded at each length. At each step, we have how many symbols remain to
31 be coded, what the last code length used was, and how many bit patterns of
32 that length remain unused. Then we add one to the code length and double the
33 number of unused patterns to graduate to the next code length. We then
34 assign all portions of the remaining symbols to that code length that
35 preserve the properties of a correct and eventually complete code. Those
36 properties are: we cannot use more bit patterns than are available; and when
37 all the symbols are used, there are exactly zero possible bit patterns
38 remaining.
39
40 The inflate Huffman decoding algorithm uses two-level lookup tables for
41 speed. There is a single first-level table to decode codes up to root bits
42 in length (root == 9 in the current inflate implementation). The table
43 has 1 << root entries and is indexed by the next root bits of input. Codes
44 shorter than root bits have replicated table entries, so that the correct
45 entry is pointed to regardless of the bits that follow the short code. If
46 the code is longer than root bits, then the table entry points to a second-
47 level table. The size of that table is determined by the longest code with
48 that root-bit prefix. If that longest code has length len, then the table
49 has size 1 << (len - root), to index the remaining bits in that set of
50 codes. Each subsequent root-bit prefix then has its own sub-table. The
51 total number of table entries required by the code is calculated
52 incrementally as the number of codes at each bit length is populated. When
53 all of the codes are shorter than root bits, then root is reduced to the
54 longest code length, resulting in a single, smaller, one-level table.
55
56 The inflate algorithm also provides for small values of root (relative to
57 the log2 of the number of symbols), where the shortest code has more bits
58 than root. In that case, root is increased to the length of the shortest
59 code. This program, by design, does not handle that case, so it is verified
60 that the number of symbols is less than 2^(root + 1).
61
62 In order to speed up the examination (by about ten orders of magnitude for
63 the default arguments), the intermediate states in the build-up of a code
64 are remembered and previously visited branches are pruned. The memory
65 required for this will increase rapidly with the total number of symbols and
66 the maximum code length in bits. However this is a very small price to pay
67 for the vast speedup.
68
69 First, all of the possible Huffman codes are counted, and reachable
70 intermediate states are noted by a non-zero count in a saved-results array.
71 Second, the intermediate states that lead to (root + 1) bit or longer codes
72 are used to look at all sub-codes from those junctures for their inflate
73 memory usage. (The amount of memory used is not affected by the number of
74 codes of root bits or less in length.) Third, the visited states in the
75 construction of those sub-codes and the associated calculation of the table
76 size is recalled in order to avoid recalculating from the same juncture.
77 Beginning the code examination at (root + 1) bit codes, which is enabled by
78 identifying the reachable nodes, accounts for about six of the orders of
79 magnitude of improvement for the default arguments. About another four
80 orders of magnitude come from not revisiting previous states. Out of
81 approximately 2x10^16 possible Huffman codes, only about 2x10^6 sub-codes
82 need to be examined to cover all of the possible table memory usage cases
83 for the default arguments of 286 symbols limited to 15-bit codes.
84
85 Note that an unsigned long long type is used for counting. It is quite easy
86 to exceed the capacity of an eight-byte integer with a large number of
87 symbols and a large maximum code length, so multiple-precision arithmetic
88 would need to replace the unsigned long long arithmetic in that case. This
89 program will abort if an overflow occurs. The big_t type identifies where
90 the counting takes place.
91
92 An unsigned long long type is also used for calculating the number of
93 possible codes remaining at the maximum length. This limits the maximum
94 code length to the number of bits in a long long minus the number of bits
95 needed to represent the symbols in a flat code. The code_t type identifies
96 where the bit pattern counting takes place.
97 */
98
99 #include <stdio.h>
100 #include <stdlib.h>
101 #include <string.h>
102 #include <assert.h>
103
104 #define local static
105
106 /* special data types */
107 typedef unsigned long long big_t; /* type for code counting */
108 typedef unsigned long long code_t; /* type for bit pattern counting */
109 struct tab { /* type for been here check */
110 size_t len; /* length of bit vector in char's */
111 char *vec; /* allocated bit vector */
112 };
113
114 /* The array for saving results, num[], is indexed with this triplet:
115
116 syms: number of symbols remaining to code
117 left: number of available bit patterns at length len
118 len: number of bits in the codes currently being assigned
119
120 Those indices are constrained thusly when saving results:
121
122 syms: 3..totsym (totsym == total symbols to code)
123 left: 2..syms - 1, but only the evens (so syms == 8 -> 2, 4, 6)
124 len: 1..max - 1 (max == maximum code length in bits)
125
126 syms == 2 is not saved since that immediately leads to a single code. left
127 must be even, since it represents the number of available bit patterns at
128 the current length, which is double the number at the previous length.
129 left ends at syms-1 since left == syms immediately results in a single code.
130 (left > sym is not allowed since that would result in an incomplete code.)
131 len is less than max, since the code completes immediately when len == max.
132
133 The offset into the array is calculated for the three indices with the
134 first one (syms) being outermost, and the last one (len) being innermost.
135 We build the array with length max-1 lists for the len index, with syms-3
136 of those for each symbol. There are totsym-2 of those, with each one
137 varying in length as a function of sym. See the calculation of index in
138 count() for the index, and the calculation of size in main() for the size
139 of the array.
140
141 For the deflate example of 286 symbols limited to 15-bit codes, the array
142 has 284,284 entries, taking up 2.17 MB for an 8-byte big_t. More than
143 half of the space allocated for saved results is actually used -- not all
144 possible triplets are reached in the generation of valid Huffman codes.
145 */
146
147 /* The array for tracking visited states, done[], is itself indexed identically
148 to the num[] array as described above for the (syms, left, len) triplet.
149 Each element in the array is further indexed by the (mem, rem) doublet,
150 where mem is the amount of inflate table space used so far, and rem is the
151 remaining unused entries in the current inflate sub-table. Each indexed
152 element is simply one bit indicating whether the state has been visited or
153 not. Since the ranges for mem and rem are not known a priori, each bit
154 vector is of a variable size, and grows as needed to accommodate the visited
155 states. mem and rem are used to calculate a single index in a triangular
156 array. Since the range of mem is expected in the default case to be about
157 ten times larger than the range of rem, the array is skewed to reduce the
158 memory usage, with eight times the range for mem than for rem. See the
159 calculations for offset and bit in beenhere() for the details.
160
161 For the deflate example of 286 symbols limited to 15-bit codes, the bit
162 vectors grow to total approximately 21 MB, in addition to the 4.3 MB done[]
163 array itself.
164 */
165
166 /* Globals to avoid propagating constants or constant pointers recursively */
167 local int max; /* maximum allowed bit length for the codes */
168 local int root; /* size of base code table in bits */
169 local int large; /* largest code table so far */
170 local size_t size; /* number of elements in num and done */
171 local int *code; /* number of symbols assigned to each bit length */
172 local big_t *num; /* saved results array for code counting */
173 local struct tab *done; /* states already evaluated array */
174
175 /* Index function for num[] and done[] */
176 #define INDEX(i,j,k) (((size_t)((i-1)>>1)*((i-2)>>1)+(j>>1)-1)*(max-1)+k-1)
177
178 /* Free allocated space. Uses globals code, num, and done. */
cleanup(void)179 local void cleanup(void)
180 {
181 size_t n;
182
183 if (done != NULL) {
184 for (n = 0; n < size; n++)
185 if (done[n].len)
186 free(done[n].vec);
187 free(done);
188 }
189 if (num != NULL)
190 free(num);
191 if (code != NULL)
192 free(code);
193 }
194
195 /* Return the number of possible Huffman codes using bit patterns of lengths
196 len through max inclusive, coding syms symbols, with left bit patterns of
197 length len unused -- return -1 if there is an overflow in the counting.
198 Keep a record of previous results in num to prevent repeating the same
199 calculation. Uses the globals max and num. */
count(int syms,int len,int left)200 local big_t count(int syms, int len, int left)
201 {
202 big_t sum; /* number of possible codes from this juncture */
203 big_t got; /* value returned from count() */
204 int least; /* least number of syms to use at this juncture */
205 int most; /* most number of syms to use at this juncture */
206 int use; /* number of bit patterns to use in next call */
207 size_t index; /* index of this case in *num */
208
209 /* see if only one possible code */
210 if (syms == left)
211 return 1;
212
213 /* note and verify the expected state */
214 assert(syms > left && left > 0 && len < max);
215
216 /* see if we've done this one already */
217 index = INDEX(syms, left, len);
218 got = num[index];
219 if (got)
220 return got; /* we have -- return the saved result */
221
222 /* we need to use at least this many bit patterns so that the code won't be
223 incomplete at the next length (more bit patterns than symbols) */
224 least = (left << 1) - syms;
225 if (least < 0)
226 least = 0;
227
228 /* we can use at most this many bit patterns, lest there not be enough
229 available for the remaining symbols at the maximum length (if there were
230 no limit to the code length, this would become: most = left - 1) */
231 most = (((code_t)left << (max - len)) - syms) /
232 (((code_t)1 << (max - len)) - 1);
233
234 /* count all possible codes from this juncture and add them up */
235 sum = 0;
236 for (use = least; use <= most; use++) {
237 got = count(syms - use, len + 1, (left - use) << 1);
238 sum += got;
239 if (got == -1 || sum < got) /* overflow */
240 return -1;
241 }
242
243 /* verify that all recursive calls are productive */
244 assert(sum != 0);
245
246 /* save the result and return it */
247 num[index] = sum;
248 return sum;
249 }
250
251 /* Return true if we've been here before, set to true if not. Set a bit in a
252 bit vector to indicate visiting this state. Each (syms,len,left) state
253 has a variable size bit vector indexed by (mem,rem). The bit vector is
254 lengthened if needed to allow setting the (mem,rem) bit. */
beenhere(int syms,int len,int left,int mem,int rem)255 local int beenhere(int syms, int len, int left, int mem, int rem)
256 {
257 size_t index; /* index for this state's bit vector */
258 size_t offset; /* offset in this state's bit vector */
259 int bit; /* mask for this state's bit */
260 size_t length; /* length of the bit vector in bytes */
261 char *vector; /* new or enlarged bit vector */
262
263 /* point to vector for (syms,left,len), bit in vector for (mem,rem) */
264 index = INDEX(syms, left, len);
265 mem -= 1 << root;
266 offset = (mem >> 3) + rem;
267 offset = ((offset * (offset + 1)) >> 1) + rem;
268 bit = 1 << (mem & 7);
269
270 /* see if we've been here */
271 length = done[index].len;
272 if (offset < length && (done[index].vec[offset] & bit) != 0)
273 return 1; /* done this! */
274
275 /* we haven't been here before -- set the bit to show we have now */
276
277 /* see if we need to lengthen the vector in order to set the bit */
278 if (length <= offset) {
279 /* if we have one already, enlarge it, zero out the appended space */
280 if (length) {
281 do {
282 length <<= 1;
283 } while (length <= offset);
284 vector = realloc(done[index].vec, length);
285 if (vector != NULL)
286 memset(vector + done[index].len, 0, length - done[index].len);
287 }
288
289 /* otherwise we need to make a new vector and zero it out */
290 else {
291 length = 1 << (len - root);
292 while (length <= offset)
293 length <<= 1;
294 vector = calloc(length, sizeof(char));
295 }
296
297 /* in either case, bail if we can't get the memory */
298 if (vector == NULL) {
299 fputs("abort: unable to allocate enough memory\n", stderr);
300 cleanup();
301 exit(1);
302 }
303
304 /* install the new vector */
305 done[index].len = length;
306 done[index].vec = vector;
307 }
308
309 /* set the bit */
310 done[index].vec[offset] |= bit;
311 return 0;
312 }
313
314 /* Examine all possible codes from the given node (syms, len, left). Compute
315 the amount of memory required to build inflate's decoding tables, where the
316 number of code structures used so far is mem, and the number remaining in
317 the current sub-table is rem. Uses the globals max, code, root, large, and
318 done. */
examine(int syms,int len,int left,int mem,int rem)319 local void examine(int syms, int len, int left, int mem, int rem)
320 {
321 int least; /* least number of syms to use at this juncture */
322 int most; /* most number of syms to use at this juncture */
323 int use; /* number of bit patterns to use in next call */
324
325 /* see if we have a complete code */
326 if (syms == left) {
327 /* set the last code entry */
328 code[len] = left;
329
330 /* complete computation of memory used by this code */
331 while (rem < left) {
332 left -= rem;
333 rem = 1 << (len - root);
334 mem += rem;
335 }
336 assert(rem == left);
337
338 /* if this is a new maximum, show the entries used and the sub-code */
339 if (mem > large) {
340 large = mem;
341 printf("max %d: ", mem);
342 for (use = root + 1; use <= max; use++)
343 if (code[use])
344 printf("%d[%d] ", code[use], use);
345 putchar('\n');
346 fflush(stdout);
347 }
348
349 /* remove entries as we drop back down in the recursion */
350 code[len] = 0;
351 return;
352 }
353
354 /* prune the tree if we can */
355 if (beenhere(syms, len, left, mem, rem))
356 return;
357
358 /* we need to use at least this many bit patterns so that the code won't be
359 incomplete at the next length (more bit patterns than symbols) */
360 least = (left << 1) - syms;
361 if (least < 0)
362 least = 0;
363
364 /* we can use at most this many bit patterns, lest there not be enough
365 available for the remaining symbols at the maximum length (if there were
366 no limit to the code length, this would become: most = left - 1) */
367 most = (((code_t)left << (max - len)) - syms) /
368 (((code_t)1 << (max - len)) - 1);
369
370 /* occupy least table spaces, creating new sub-tables as needed */
371 use = least;
372 while (rem < use) {
373 use -= rem;
374 rem = 1 << (len - root);
375 mem += rem;
376 }
377 rem -= use;
378
379 /* examine codes from here, updating table space as we go */
380 for (use = least; use <= most; use++) {
381 code[len] = use;
382 examine(syms - use, len + 1, (left - use) << 1,
383 mem + (rem ? 1 << (len - root) : 0), rem << 1);
384 if (rem == 0) {
385 rem = 1 << (len - root);
386 mem += rem;
387 }
388 rem--;
389 }
390
391 /* remove entries as we drop back down in the recursion */
392 code[len] = 0;
393 }
394
395 /* Look at all sub-codes starting with root + 1 bits. Look at only the valid
396 intermediate code states (syms, left, len). For each completed code,
397 calculate the amount of memory required by inflate to build the decoding
398 tables. Find the maximum amount of memory required and show the code that
399 requires that maximum. Uses the globals max, root, and num. */
enough(int syms)400 local void enough(int syms)
401 {
402 int n; /* number of remaing symbols for this node */
403 int left; /* number of unused bit patterns at this length */
404 size_t index; /* index of this case in *num */
405
406 /* clear code */
407 for (n = 0; n <= max; n++)
408 code[n] = 0;
409
410 /* look at all (root + 1) bit and longer codes */
411 large = 1 << root; /* base table */
412 if (root < max) /* otherwise, there's only a base table */
413 for (n = 3; n <= syms; n++)
414 for (left = 2; left < n; left += 2)
415 {
416 /* look at all reachable (root + 1) bit nodes, and the
417 resulting codes (complete at root + 2 or more) */
418 index = INDEX(n, left, root + 1);
419 if (root + 1 < max && num[index]) /* reachable node */
420 examine(n, root + 1, left, 1 << root, 0);
421
422 /* also look at root bit codes with completions at root + 1
423 bits (not saved in num, since complete), just in case */
424 if (num[index - 1] && n <= left << 1)
425 examine((n - left) << 1, root + 1, (n - left) << 1,
426 1 << root, 0);
427 }
428
429 /* done */
430 printf("done: maximum of %d table entries\n", large);
431 }
432
433 /*
434 Examine and show the total number of possible Huffman codes for a given
435 maximum number of symbols, initial root table size, and maximum code length
436 in bits -- those are the command arguments in that order. The default
437 values are 286, 9, and 15 respectively, for the deflate literal/length code.
438 The possible codes are counted for each number of coded symbols from two to
439 the maximum. The counts for each of those and the total number of codes are
440 shown. The maximum number of inflate table entires is then calculated
441 across all possible codes. Each new maximum number of table entries and the
442 associated sub-code (starting at root + 1 == 10 bits) is shown.
443
444 To count and examine Huffman codes that are not length-limited, provide a
445 maximum length equal to the number of symbols minus one.
446
447 For the deflate literal/length code, use "enough". For the deflate distance
448 code, use "enough 30 6".
449
450 This uses the %llu printf format to print big_t numbers, which assumes that
451 big_t is an unsigned long long. If the big_t type is changed (for example
452 to a multiple precision type), the method of printing will also need to be
453 updated.
454 */
main(int argc,char ** argv)455 int main(int argc, char **argv)
456 {
457 int syms; /* total number of symbols to code */
458 int n; /* number of symbols to code for this run */
459 big_t got; /* return value of count() */
460 big_t sum; /* accumulated number of codes over n */
461
462 /* set up globals for cleanup() */
463 code = NULL;
464 num = NULL;
465 done = NULL;
466
467 /* get arguments -- default to the deflate literal/length code */
468 syms = 286;
469 root = 9;
470 max = 15;
471 if (argc > 1) {
472 syms = atoi(argv[1]);
473 if (argc > 2) {
474 root = atoi(argv[2]);
475 if (argc > 3)
476 max = atoi(argv[3]);
477 }
478 }
479 if (argc > 4 || syms < 2 || root < 1 || max < 1) {
480 fputs("invalid arguments, need: [sym >= 2 [root >= 1 [max >= 1]]]\n",
481 stderr);
482 return 1;
483 }
484
485 /* if not restricting the code length, the longest is syms - 1 */
486 if (max > syms - 1)
487 max = syms - 1;
488
489 /* determine the number of bits in a code_t */
490 n = 0;
491 while (((code_t)1 << n) != 0)
492 n++;
493
494 /* make sure that the calculation of most will not overflow */
495 if (max > n || syms - 2 >= (((code_t)0 - 1) >> (max - 1))) {
496 fputs("abort: code length too long for internal types\n", stderr);
497 return 1;
498 }
499
500 /* reject impossible code requests */
501 if (syms - 1 > ((code_t)1 << max) - 1) {
502 fprintf(stderr, "%d symbols cannot be coded in %d bits\n",
503 syms, max);
504 return 1;
505 }
506
507 /* allocate code vector */
508 code = calloc(max + 1, sizeof(int));
509 if (code == NULL) {
510 fputs("abort: unable to allocate enough memory\n", stderr);
511 return 1;
512 }
513
514 /* determine size of saved results array, checking for overflows,
515 allocate and clear the array (set all to zero with calloc()) */
516 if (syms == 2) /* iff max == 1 */
517 num = NULL; /* won't be saving any results */
518 else {
519 size = syms >> 1;
520 if (size > ((size_t)0 - 1) / (n = (syms - 1) >> 1) ||
521 (size *= n, size > ((size_t)0 - 1) / (n = max - 1)) ||
522 (size *= n, size > ((size_t)0 - 1) / sizeof(big_t)) ||
523 (num = calloc(size, sizeof(big_t))) == NULL) {
524 fputs("abort: unable to allocate enough memory\n", stderr);
525 cleanup();
526 return 1;
527 }
528 }
529
530 /* count possible codes for all numbers of symbols, add up counts */
531 sum = 0;
532 for (n = 2; n <= syms; n++) {
533 got = count(n, 1, 2);
534 sum += got;
535 if (got == -1 || sum < got) { /* overflow */
536 fputs("abort: can't count that high!\n", stderr);
537 cleanup();
538 return 1;
539 }
540 printf("%llu %d-codes\n", got, n);
541 }
542 printf("%llu total codes for 2 to %d symbols", sum, syms);
543 if (max < syms - 1)
544 printf(" (%d-bit length limit)\n", max);
545 else
546 puts(" (no length limit)");
547
548 /* allocate and clear done array for beenhere() */
549 if (syms == 2)
550 done = NULL;
551 else if (size > ((size_t)0 - 1) / sizeof(struct tab) ||
552 (done = calloc(size, sizeof(struct tab))) == NULL) {
553 fputs("abort: unable to allocate enough memory\n", stderr);
554 cleanup();
555 return 1;
556 }
557
558 /* find and show maximum inflate table usage */
559 if (root > max) /* reduce root to max length */
560 root = max;
561 if (syms < ((code_t)1 << (root + 1)))
562 enough(syms);
563 else
564 puts("cannot handle minimum code lengths > root");
565
566 /* done */
567 cleanup();
568 return 0;
569 }
570