5e0fa7faff227beb16301f80aac3c95232f0b9e4
[wimlib] / src / lz_sarray.c
1 /*
2  * lz_sarray.c
3  *
4  * Suffix array match-finder for Lempel-Ziv compression.
5  */
6
7 /*
8  * Copyright (c) 2013, 2014 Eric Biggers.  All rights reserved.
9  *
10  * Redistribution and use in source and binary forms, with or without
11  * modification, are permitted provided that the following conditions
12  * are met:
13  *
14  * 1. Redistributions of source code must retain the above copyright
15  *    notice, this list of conditions and the following disclaimer.
16  *
17  * 2. Redistributions in binary form must reproduce the above copyright
18  *    notice, this list of conditions and the following disclaimer in the
19  *    documentation and/or other materials provided with the distribution.
20  *
21  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS "AS IS" AND
22  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
23  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
24  * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE
25  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
26  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
27  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
28  * BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
29  * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
30  * OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
31  * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32  */
33
34 #ifdef HAVE_CONFIG_H
35 #  include "config.h"
36 #endif
37
38 #include "wimlib/lz_sarray.h"
39 #include "wimlib/util.h"
40 #include "divsufsort/divsufsort.h"
41 #include <string.h>
42
43 /* If ENABLE_LZ_DEBUG is defined, verify that the suffix array satisfies its
44  * definition.
45  *
46  * @SA          The constructed suffix array.
47  * @T           The original data.
48  * @found       Temporary 'bool' array of length @n.
49  * @n           Length of the data (length of @SA, @T, and @found arrays).
50  *
51  * WARNING: this is for debug use only as it does not necessarily run in linear
52  * time!!!  */
53 static void
54 verify_suffix_array(const lz_sarray_pos_t SA[restrict],
55                     const u8 T[restrict],
56                     bool found[restrict],
57                     lz_sarray_pos_t n)
58 {
59 #ifdef ENABLE_LZ_DEBUG
60         /* Ensure the SA contains exactly one of each i in [0, n - 1].  */
61         for (lz_sarray_pos_t i = 0; i < n; i++)
62                 found[i] = false;
63         for (lz_sarray_pos_t r = 0; r < n; r++) {
64                 lz_sarray_pos_t i = SA[r];
65                 LZ_ASSERT(i < n);
66                 LZ_ASSERT(!found[i]);
67                 found[i] = true;
68         }
69
70         /* Ensure the suffix with rank r is lexicographically lesser than the
71          * suffix with rank (r + 1) for all r in [0, n - 2].  */
72         for (lz_sarray_pos_t r = 0; r < n - 1; r++) {
73
74                 lz_sarray_pos_t i1 = SA[r];
75                 lz_sarray_pos_t i2 = SA[r + 1];
76
77                 lz_sarray_pos_t n1 = n - i1;
78                 lz_sarray_pos_t n2 = n - i2;
79
80                 int res = memcmp(&T[i1], &T[i2], min(n1, n2));
81                 LZ_ASSERT(res < 0 || (res == 0 && n1 < n2));
82         }
83 #endif /* ENABLE_LZ_DEBUG  */
84 }
85
86 /* Compute the inverse suffix array @ISA from the suffix array @SA in linear
87  * time.
88  *
89  * Whereas the suffix array is a mapping from suffix rank to suffix position,
90  * the inverse suffix array is a mapping from suffix position to suffix rank.
91  */
92 static void
93 compute_inverse_suffix_array(lz_sarray_pos_t ISA[restrict],
94                              const lz_sarray_pos_t SA[restrict],
95                              lz_sarray_pos_t n)
96 {
97         lz_sarray_pos_t r;
98
99         for (r = 0; r < n; r++)
100                 ISA[SA[r]] = r;
101 }
102
103
104 /* Compute the LCP (Longest Common Prefix) array in linear time.
105  *
106  * LCP[r] will be the length of the longest common prefix between the suffixes
107  * with positions SA[r - 1] and  SA[r].  LCP[0] will be undefined.
108  *
109  * Algorithm adapted from Kasai et al. 2001: "Linear-Time Longest-Common-Prefix
110  * Computation in Suffix Arrays and Its Applications".  Modified slightly to
111  * take into account that with bytes in the real world, there is no unique
112  * symbol at the end of the string.  */
113 static void
114 compute_lcp_array(lz_sarray_pos_t LCP[restrict],
115                   const lz_sarray_pos_t SA[restrict],
116                   const lz_sarray_pos_t ISA[restrict],
117                   const u8 T[restrict],
118                   lz_sarray_pos_t n)
119 {
120         lz_sarray_pos_t h, i, r, j, lim;
121
122         h = 0;
123         for (i = 0; i < n; i++) {
124                 r = ISA[i];
125                 if (r > 0) {
126                         j = SA[r - 1];
127                         lim = min(n - i, n - j);
128
129                         while (h < lim && T[i + h] == T[j + h])
130                                 h++;
131                         LCP[r] = h;
132                         if (h > 0)
133                                 h--;
134                 }
135         }
136 }
137
138 /* If ENABLE_LZ_DEBUG is defined, verify that the LCP (Longest Common Prefix)
139  * array satisfies its definition.
140  *
141  * WARNING: this is for debug use only as it does not necessarily run in linear
142  * time!!!  */
143 static void
144 verify_lcp_array(lz_sarray_pos_t LCP[restrict],
145                  const lz_sarray_pos_t SA[restrict],
146                  const u8 T[restrict],
147                  lz_sarray_pos_t n)
148 {
149 #ifdef ENABLE_LZ_DEBUG
150         for (lz_sarray_pos_t r = 0; r < n - 1; r++) {
151                 lz_sarray_pos_t i1 = SA[r];
152                 lz_sarray_pos_t i2 = SA[r + 1];
153                 lz_sarray_pos_t lcp = LCP[r + 1];
154
155                 lz_sarray_pos_t n1 = n - i1;
156                 lz_sarray_pos_t n2 = n - i2;
157
158                 LZ_ASSERT(lcp <= min(n1, n2));
159
160                 LZ_ASSERT(memcmp(&T[i1], &T[i2], lcp) == 0);
161                 if (lcp < min(n1, n2))
162                         LZ_ASSERT(T[i1 + lcp] != T[i2 + lcp]);
163         }
164 #endif /* ENABLE_LZ_DEBUG */
165 }
166
167 /* Initialize the SA link array in linear time.
168  *
169  * This is similar to computing the LPF (Longest Previous Factor) array, which
170  * is addressed in several papers.  In particular the algorithms below are based
171  * on Crochemore et al. 2009: "LPF computation revisited".  However, this
172  * match-finder does not actually compute or use the LPF array per se.  Rather,
173  * this function sets up some information necessary to compute the LPF array,
174  * but later lz_sarray_get_matches() actually uses this information to search
175  * the suffix array directly and can keep searching beyond the first (longest)
176  * match whose length would be placed in the LPF array.  This difference from
177  * the theoretical work is necessary because in many real compression formats
178  * matches take variable numbers of bits to encode, so a decent parser needs to
179  * consider more than just the longest match with unspecified offset.
180  *
181  * Note: We cap the lcpprev and lcpnext values to the maximum match length so
182  * that the match-finder need not worry about it later, in the inner loop.
183  */
184 static void
185 init_salink(struct salink link[restrict],
186             const lz_sarray_pos_t LCP[restrict],
187             const lz_sarray_pos_t SA[restrict],
188             const u8 T[restrict],
189             lz_sarray_pos_t n,
190             lz_sarray_len_t max_match_len)
191 {
192         /* Compute salink.next and salink.lcpnext.  */
193         link[n - 1].next = LZ_SARRAY_POS_MAX;
194         link[n - 1].lcpnext = 0;
195         for (lz_sarray_pos_t r = n - 2; r != LZ_SARRAY_POS_MAX; r--) {
196                 lz_sarray_pos_t t = r + 1;
197                 lz_sarray_pos_t l = LCP[t];
198                 while (t != LZ_SARRAY_POS_MAX && SA[t] > SA[r]) {
199                         l = min(l, link[t].lcpnext);
200                         t = link[t].next;
201                 }
202                 link[r].next = t;
203                 link[r].lcpnext = min(l, max_match_len);
204         }
205
206         /* Compute salink.prev and salink.lcpprev.  */
207         link[0].prev = LZ_SARRAY_POS_MAX;
208         link[0].lcpprev = 0;
209         for (lz_sarray_pos_t r = 1; r < n; r++) {
210                 lz_sarray_pos_t t = r - 1;
211                 lz_sarray_pos_t l = LCP[r];
212                 while (t != LZ_SARRAY_POS_MAX && SA[t] > SA[r]) {
213                         l = min(l, link[t].lcpprev);
214                         t = link[t].prev;
215                 }
216                 link[r].prev = t;
217                 link[r].lcpprev = min(l, max_match_len);
218         }
219 }
220
221 /* If ENABLE_LZ_DEBUG is defined, verify the values computed by init_salink().
222  *
223  * WARNING: this is for debug use only as it does not necessarily run in linear
224  * time!!!  */
225 static void
226 verify_salink(const struct salink link[],
227               const lz_sarray_pos_t SA[],
228               const u8 T[],
229               lz_sarray_pos_t n,
230               lz_sarray_len_t max_match_len)
231 {
232 #ifdef ENABLE_LZ_DEBUG
233         for (lz_sarray_pos_t r = 0; r < n; r++) {
234                 for (lz_sarray_pos_t prev = r; ; ) {
235                         if (prev == 0) {
236                                 LZ_ASSERT(link[r].prev == LZ_SARRAY_POS_MAX);
237                                 LZ_ASSERT(link[r].lcpprev == 0);
238                                 break;
239                         }
240
241                         prev--;
242
243                         if (SA[prev] < SA[r]) {
244                                 LZ_ASSERT(link[r].prev == prev);
245                                 LZ_ASSERT(link[r].lcpprev <= n - SA[prev]);
246                                 LZ_ASSERT(link[r].lcpprev <= n - SA[r]);
247                                 LZ_ASSERT(link[r].lcpprev <= max_match_len);
248                                 LZ_ASSERT(0 == memcmp(&T[SA[prev]],
249                                                       &T[SA[r]],
250                                                       link[r].lcpprev));
251                                 if (link[r].lcpprev < n - SA[prev] &&
252                                     link[r].lcpprev < n - SA[r] &&
253                                     link[r].lcpprev < max_match_len)
254                                 {
255                                         LZ_ASSERT(T[SA[prev] + link[r].lcpprev] !=
256                                                   T[SA[r] + link[r].lcpprev]);
257                                 }
258                                 break;
259                         }
260                 }
261
262                 for (lz_sarray_pos_t next = r; ; ) {
263                         if (next == n - 1) {
264                                 LZ_ASSERT(link[r].next == LZ_SARRAY_POS_MAX);
265                                 LZ_ASSERT(link[r].lcpnext == 0);
266                                 break;
267                         }
268
269                         next++;
270
271                         if (SA[next] < SA[r]) {
272                                 LZ_ASSERT(link[r].next == next);
273                                 LZ_ASSERT(link[r].lcpnext <= n - SA[next]);
274                                 LZ_ASSERT(link[r].lcpnext <= n - SA[r]);
275                                 LZ_ASSERT(link[r].lcpnext <= max_match_len);
276                                 LZ_ASSERT(0 == memcmp(&T[SA[next]],
277                                                       &T[SA[r]],
278                                                       link[r].lcpnext));
279                                 if (link[r].lcpnext < n - SA[next] &&
280                                     link[r].lcpnext < n - SA[r] &&
281                                     link[r].lcpnext < max_match_len)
282                                 {
283                                         LZ_ASSERT(T[SA[next] + link[r].lcpnext] !=
284                                                   T[SA[r] + link[r].lcpnext]);
285
286                                 }
287                                 break;
288                         }
289                 }
290         }
291 #endif
292 }
293
294 /*
295  * Initialize the suffix array match-finder.
296  *
297  * @mf
298  *      The suffix array match-finder structure to initialize.  This structure
299  *      is expected to be zeroed before this function is called.  In the case
300  *      that this function fails, lz_sarray_destroy() should be called to free
301  *      any memory that may have been allocated.
302  *
303  * @max_window_size
304  *      The maximum window size to support.  This must be greater than 0.
305  *
306  *      The amount of needed memory will depend on this value; see
307  *      lz_sarray_get_needed_memory() for details.
308  *
309  * @min_match_len
310  *      The minimum length of each match to be found.  Must be greater than 0.
311  *
312  * @max_match_len
313  *      The maximum length of each match to be found.  Must be greater than or
314  *      equal to @min_match_len.
315  *
316  * @max_matches_to_consider
317  *      The maximum number of matches to consider at each position.  This should
318  *      be greater than @max_matches_to_return because @max_matches_to_consider
319  *      counts all the returned matches as well as matches of equal length to
320  *      returned matches that were not returned.  This parameter bounds the
321  *      amount of work the match-finder does at any one position.  This could be
322  *      anywhere from 1 to 100+ depending on the compression ratio and
323  *      performance desired.
324  *
325  * @max_matches_to_return
326  *      Maximum number of matches to return at each position.  Because of the
327  *      suffix array search algorithm, the order in which matches are returned
328  *      will be from longest to shortest, so cut-offs due to this parameter will
329  *      only result in shorter matches being discarded.  This parameter could be
330  *      anywhere from 1 to (@max_match_len - @min_match_len + 1) depending on
331  *      the compression performance desired.  However, making it even moderately
332  *      large (say, greater than 3) may not be very helpful due to the property
333  *      that the matches are returned from longest to shortest.  But the main
334  *      thing to keep in mind is that if the compressor decides to output a
335  *      shorter-than-possible match, ideally it would be best to choose the best
336  *      match of the desired length rather than truncate a longer match to that
337  *      length.
338  *
339  * After initialization, the suffix-array match-finder can be used for any
340  * number of input strings (windows) of length less than or equal to
341  * @max_window_size by successive calls to lz_sarray_load_window().
342  *
343  * Returns %true on success, or %false if sufficient memory could not be
344  * allocated.  See the note for @max_window_size above regarding the needed
345  * memory size.
346  */
347 bool
348 lz_sarray_init(struct lz_sarray *mf,
349                lz_sarray_pos_t max_window_size,
350                lz_sarray_len_t min_match_len,
351                lz_sarray_len_t max_match_len,
352                u32 max_matches_to_consider,
353                u32 max_matches_to_return)
354 {
355         LZ_ASSERT(min_match_len > 0);
356         LZ_ASSERT(max_window_size > 0);
357         LZ_ASSERT(max_match_len >= min_match_len);
358
359         mf->max_window_size = max_window_size;
360         mf->min_match_len = min_match_len;
361         mf->max_match_len = max_match_len;
362         mf->max_matches_to_consider = max_matches_to_consider;
363         mf->max_matches_to_return = max_matches_to_return;
364
365         /* SA and ISA will share the same storage block.  */
366         if ((u64)2 * max_window_size * sizeof(mf->SA[0]) !=
367                  2 * max_window_size * sizeof(mf->SA[0]))
368                 return false;
369         mf->SA = MALLOC(2 * max_window_size * sizeof(mf->SA[0]));
370         if (mf->SA == NULL)
371                 return false;
372
373         if ((u64)max_window_size * sizeof(mf->salink[0]) !=
374                  max_window_size * sizeof(mf->salink[0]))
375                 return false;
376         mf->salink = MALLOC(max_window_size * sizeof(mf->salink[0]));
377         if (mf->salink == NULL)
378                 return false;
379
380         return true;
381 }
382
383 /*
384  * Return the number of bytes of memory that lz_sarray_init() would allocate for
385  * the specified maximum window size.
386  *
387  * This should be (20 * @max_window_size) unless the type definitions have been
388  * changed.
389  */
390 u64
391 lz_sarray_get_needed_memory(lz_sarray_pos_t max_window_size)
392 {
393         u64 size = 0;
394
395         /* SA and ISA: 8 bytes per position  */
396         size += (u64)max_window_size * 2 * sizeof(((struct lz_sarray*)0)->SA[0]);
397
398         /* salink: 12 bytes per position  */
399         size += (u64)max_window_size * sizeof(((struct lz_sarray*)0)->salink[0]);
400
401         return size;
402 }
403
404 /*
405  * Prepare the suffix array match-finder to scan the specified window for
406  * matches.
407  *
408  * @mf  Suffix array match-finder previously initialized with lz_sarray_init().
409  *
410  * @T   Window, or "block", in which to find matches.
411  *
412  * @n   Size of window in bytes.  This must be positive and less than or equal
413  *      to the @max_window_size passed to lz_sarray_init().
414  *
415  * This function runs in linear time (relative to @n).
416  */
417 void
418 lz_sarray_load_window(struct lz_sarray *mf, const u8 T[], lz_sarray_pos_t n)
419 {
420         lz_sarray_pos_t *ISA, *LCP;
421
422         LZ_ASSERT(n > 0 && n <= mf->max_window_size);
423
424         /* Compute SA (Suffix Array).
425          *
426          * divsufsort() needs temporary space --- one array with 256 spaces and
427          * one array with 65536 spaces.  The implementation has been modified
428          * from the original to use the provided temporary space instead of
429          * allocating its own.
430          *
431          * We also check at build-time that divsufsort() uses the same integer
432          * size expected by this code.  Unfortunately, divsufsort breaks if
433          * 'sa_idx_t' is defined to be a 16-bit integer; however, that would
434          * limit blocks to only 65536 bytes anyway.  */
435         LZ_ASSERT(mf->max_window_size * sizeof(mf->SA[0])
436                   >= 256 * sizeof(saidx_t));
437         LZ_ASSERT(mf->max_window_size * sizeof(mf->salink[0])
438                   >= 256 * 256 * sizeof(saidx_t));
439         BUILD_BUG_ON(sizeof(lz_sarray_pos_t) != sizeof(saidx_t));
440
441         divsufsort(T, mf->SA, n, (saidx_t*)&mf->SA[n], (saidx_t*)mf->salink);
442
443         BUILD_BUG_ON(sizeof(bool) > sizeof(mf->salink[0]));
444         verify_suffix_array(mf->SA, T, (bool*)mf->salink, n);
445
446         /* Compute ISA (Inverse Suffix Array) in a preliminary position.
447          *
448          * This is just a trick to save memory.  Since LCP is unneeded after
449          * this function, it can be computed in any available space.  The
450          * storage for the ISA is the best choice because the ISA can be built
451          * quickly in salink for now, then re-built in its real location at the
452          * end.  This is probably worth it because computing the ISA from the SA
453          * is very fast, and since this match-finder is memory-hungry we'd like
454          * to save as much memory as possible.  */
455         ISA = (lz_sarray_pos_t*)mf->salink;
456         compute_inverse_suffix_array(ISA, mf->SA, n);
457
458         /* Compute LCP (Longest Common Prefix) array.  */
459         LCP = mf->SA + n;
460         compute_lcp_array(LCP, mf->SA, ISA, T, n);
461         verify_lcp_array(LCP, mf->SA, T, n);
462
463         /* Initialize suffix array links.  */
464         init_salink(mf->salink, LCP, mf->SA, T, n, mf->max_match_len);
465         verify_salink(mf->salink, mf->SA, T, n, mf->max_match_len);
466
467         /* Compute ISA (Inverse Suffix Array) in its final position.  */
468         ISA = mf->SA + n;
469         compute_inverse_suffix_array(ISA, mf->SA, n);
470
471         /* Save new variables and return.  */
472         mf->ISA = ISA;
473         mf->cur_pos = 0;
474         mf->window_size = n;
475 }
476
477 /* Free memory allocated for the suffix array match-finder.  */
478 void
479 lz_sarray_destroy(struct lz_sarray *mf)
480 {
481         FREE(mf->SA);
482         FREE(mf->salink);
483 }