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