Compare commits

...
Author SHA1 Message Date
ddidderr 4f6a8f0185 feat(rust): port divsufsort
Move the dictionary builder's suffix-array construction from
lib/dictBuilder/divsufsort.c to rust/src/divsufsort.rs, the first
dictBuilder module to migrate.  It rides on the dict-builder cargo
feature dimension introduced by the previous commit.

divsufsort() is a self-contained algorithm (two-stage sort of type-B*
substrings via sssort, rank refinement via trsort, then induced sorting
of the full array), so its context-free signature allows a direct symbol
takeover: the Rust #[no_mangle] export provides the existing `divsufsort`
symbol and the C file becomes a declaration-only shim that just keeps the
header's prototypes in the build.  Only divsufsort() moved; divbwt() has
no callers anywhere in zstd, so it is now declaration-only, keeping the
Rust export surface minimal.  The unused openMP parameter is retained for
signature compatibility (zstd never defines LIBBSC_OPENMP).

The port is a mechanical translation of the exact configuration zstd
compiles: ALPHABET_SIZE=256, SS_INSERTIONSORT_THRESHOLD=8,
SS_BLOCKSIZE=1024, SS_MISORT_STACKSIZE=16, SS_SMERGE_STACKSIZE=32,
TR_STACKSIZE=64.  Every C `int*` cursor into the SA buffer becomes an
`isize` index into a single `&mut [i32]` slice, preserving the pointer
arithmetic (including transient one-before-the-range cursors and the
bitwise-complement rank marking) while staying bounds-checked; all value
arithmetic keeps C int semantics.  The C -1/-2 error results are
preserved, with Vec::try_reserve_exact standing in for the bucket-array
malloc failure path.  Behavior is bit-identical by construction and by
measurement (see test plan); runtime on an 11 MB training buffer is
within ~5% of the C build end-to-end.

Users see no behavioral change: dictionaries trained through
ZDICT_trainFromBuffer_legacy() are byte-identical to the C build.  The
only external difference is that the never-called `divbwt` symbol is no
longer defined in the library.

Test plan:
- cd rust && cargo fmt --check && cargo clippy --all-targets
  -- -D warnings && cargo test --all-targets && cargo build --release
  (125 tests pass; new unit tests cover empty/one/two-byte inputs,
  all-equal bytes, an exact hand-computed "abracadabra" SA, and
  fixed-seed LCG buffers at 256-, 4-, and 2-symbol alphabets verified
  against a naive reference sort plus permutation/sorted invariants)
- Feature matrix: cargo build --release --no-default-features
  --features compression,decompression (and decompression-only,
  compression-only, compression,dict-builder); `divsufsort` is exported
  only when dict-builder is enabled
- make -C tests fuzzer && ./tests/fuzzer -i1 --no-big-tests (includes
  ZDICT training tests): pass
- make -C tests test-rust-lib-smoke: pass
- make -C tests test-invalidDictionaries: pass
- make -C programs zstd zstd-dictBuilder zstd-small zstd-compress
  zstd-decompress: build; compress/decompress round-trip verified
- Byte-identity vs pristine C build (commit 959e4852): a harness calling
  ZDICT_trainFromBuffer_legacy() (the only zstd path reaching
  divsufsort) and divsufsort() directly, linked against both libzstd.a
  builds, produces byte-identical dictionaries (80,288 B and full
  112,640 B capacity) and byte-identical suffix arrays on a 1 MB source
  set and an 11 MB binary/repetitive set; a differential driver over 148
  random and structured buffers (sizes 3..6000, alphabets 1..256,
  Fibonacci word, sawtooth, 6 KB near-constant) shows zero mismatches.
  The CLI --train path could not be exercised because the Rust CLI
  frontend rejects --train in both the pristine and ported builds (a
  pre-existing migration gap unrelated to this change).
2026-07-11 14:21:09 +02:00
ddidderr eecd2e168a build(rust): add dict-builder cargo feature dimension
The dictionary-builder sources (lib/dictBuilder) are about to start moving
to Rust, beginning with divsufsort. The Rust crate previously only modeled
the compression/decompression module split plus the forced-HUF decoder
modes, so no build could express "this C configuration includes (or
excludes) dictBuilder" to Cargo. Without that, a Rust archive could carry
dictBuilder modules into a build whose C side disabled them, or worse,
omit a migrated implementation from a build whose C shims require it.

Add a `dict-builder` cargo feature and thread it through every build that
consumes the Rust static archive, mirroring exactly how each build system
already gates the dictBuilder C sources:

- rust/Cargo.toml: new `dict-builder` feature, included in the default
  set because the C library builds dictBuilder by default
  (ZSTD_LIB_DICTBUILDER ?= 1). The feature is empty until the first
  dictBuilder module lands.
- lib/Makefile: RUST_CARGO_FEATURES gains dict-builder when
  ZSTD_LIB_DICTBUILDER is enabled, following the existing
  ZSTD_LIB_COMPRESSION/ZSTD_LIB_DECOMPRESSION pattern. The archive
  directory naming grows a matching `b<0|1>` dimension
  (c1-d1-b1-default etc.) so differently configured archives never
  collide; the repeated config prefix is factored into
  RUST_MODULE_CONFIG.
- programs/Makefile: the full-featured archives now request
  compression,decompression,dict-builder (equal to the default set, so
  the target directory stays shared with tests). The partial-library
  variants gain the `b0` name dimension, and zstd-dictBuilder gets its
  own lib-c1-d0-b1 archive because it compiles the dictBuilder C sources
  without decompression; it previously shared the compression-only
  archive, which will lack the migrated dictBuilder symbols.
- tests/Makefile: no flag change needed since tests use the crate default
  feature set; a comment now records that dict-builder arrives that way.
- build/cmake/lib/CMakeLists.txt: ZSTD_BUILD_DICTBUILDER now adds the
  dict-builder feature and a `b<0|1>` component in the Rust build-config
  directory name, in lockstep with the DictBuilderSources gating.
- build/meson/lib/meson.build: meson compiles the dictBuilder sources
  unconditionally, so the feature list and config name gain dict-builder
  unconditionally (c1-d1-b1-<huf-mode>).

The `dict-builder` feature deliberately does not imply `compression`.
lib/Makefile forces ZSTD_LIB_DICTBUILDER=0 when compression is disabled,
but CMake does not couple the two options, so encoding the C-side
constraint in Cargo would make the Rust archive diverge from the C source
list in that (already unsupported) CMake configuration.

Test plan:
- cd rust && cargo build --release
- cargo build --release --no-default-features \
    --features compression,decompression
- cargo build --release --no-default-features \
    --features compression,dict-builder
- Full validation (fuzzer, smoke tests, dictionary byte-identity) runs
  with the follow-up commit that ports divsufsort onto this scaffolding.
2026-07-11 09:38:04 +02:00
10 changed files with 2811 additions and 1903 deletions
+6 -1
View File
@@ -142,6 +142,7 @@ endif()
set(_zstd_rust_features) set(_zstd_rust_features)
set(_zstd_rust_compression 0) set(_zstd_rust_compression 0)
set(_zstd_rust_decompression 0) set(_zstd_rust_decompression 0)
set(_zstd_rust_dictbuilder 0)
if(ZSTD_BUILD_COMPRESSION) if(ZSTD_BUILD_COMPRESSION)
list(APPEND _zstd_rust_features compression) list(APPEND _zstd_rust_features compression)
set(_zstd_rust_compression 1) set(_zstd_rust_compression 1)
@@ -150,6 +151,10 @@ if(ZSTD_BUILD_DECOMPRESSION)
list(APPEND _zstd_rust_features decompression) list(APPEND _zstd_rust_features decompression)
set(_zstd_rust_decompression 1) set(_zstd_rust_decompression 1)
endif() endif()
if(ZSTD_BUILD_DICTBUILDER)
list(APPEND _zstd_rust_features dict-builder)
set(_zstd_rust_dictbuilder 1)
endif()
set(_zstd_rust_huf_mode default) set(_zstd_rust_huf_mode default)
if(_zstd_huf_force_x1) if(_zstd_huf_force_x1)
@@ -190,7 +195,7 @@ if(_zstd_rust_features)
endif() endif()
set(_zstd_rust_build_config set(_zstd_rust_build_config
"c${_zstd_rust_compression}-d${_zstd_rust_decompression}-${_zstd_rust_huf_mode}") "c${_zstd_rust_compression}-d${_zstd_rust_decompression}-b${_zstd_rust_dictbuilder}-${_zstd_rust_huf_mode}")
set(ZSTD_RUST_MANIFEST "${ZSTD_SOURCE_DIR}/rust/Cargo.toml") set(ZSTD_RUST_MANIFEST "${ZSTD_SOURCE_DIR}/rust/Cargo.toml")
set(ZSTD_RUST_TARGET_DIR set(ZSTD_RUST_TARGET_DIR
"${CMAKE_CURRENT_BINARY_DIR}/rust-target/${_zstd_rust_build_config}") "${CMAKE_CURRENT_BINARY_DIR}/rust-target/${_zstd_rust_build_config}")
+4 -2
View File
@@ -73,7 +73,9 @@ if rust_huf_force_x1 and rust_huf_force_x2
error('HUF_FORCE_DECOMPRESS_X1 and HUF_FORCE_DECOMPRESS_X2 are mutually exclusive') error('HUF_FORCE_DECOMPRESS_X1 and HUF_FORCE_DECOMPRESS_X2 are mutually exclusive')
endif endif
rust_features = ['compression', 'decompression'] # Meson always compiles the dictBuilder sources above, so the Rust archive
# must always carry the matching dict-builder module set.
rust_features = ['compression', 'decompression', 'dict-builder']
rust_huf_mode = 'default' rust_huf_mode = 'default'
rust_huf_c_args = [] rust_huf_c_args = []
if rust_huf_force_x1 if rust_huf_force_x1
@@ -96,7 +98,7 @@ if rust_target == ''
endif endif
endif endif
rust_build_config = 'c1-d1-' + rust_huf_mode rust_build_config = 'c1-d1-b1-' + rust_huf_mode
rust_target_dir = join_paths(meson.current_build_dir(), 'rust-target', rust_build_config) rust_target_dir = join_paths(meson.current_build_dir(), 'rust-target', rust_build_config)
is_msvc = cc_id == compiler_msvc or cc_id == 'clang-cl' is_msvc = cc_id == compiler_msvc or cc_id == 'clang-cl'
rust_staticlib_name = is_msvc ? 'zstd_rs.lib' : 'libzstd_rs.a' rust_staticlib_name = is_msvc ? 'zstd_rs.lib' : 'libzstd_rs.a'
+7 -3
View File
@@ -84,16 +84,20 @@ endif
ifneq ($(ZSTD_LIB_DECOMPRESSION),0) ifneq ($(ZSTD_LIB_DECOMPRESSION),0)
RUST_CARGO_FEATURES += decompression RUST_CARGO_FEATURES += decompression
endif endif
ifneq ($(ZSTD_LIB_DICTBUILDER),0)
RUST_CARGO_FEATURES += dict-builder
endif
RUST_MODULE_CONFIG := c$(ZSTD_LIB_COMPRESSION)-d$(ZSTD_LIB_DECOMPRESSION)-b$(ZSTD_LIB_DICTBUILDER)
RUST_HUF_FEATURE := RUST_HUF_FEATURE :=
RUST_BUILD_CONFIG := c$(ZSTD_LIB_COMPRESSION)-d$(ZSTD_LIB_DECOMPRESSION)-default RUST_BUILD_CONFIG := $(RUST_MODULE_CONFIG)-default
ifneq ($(RUST_HUF_FORCE_X1),0) ifneq ($(RUST_HUF_FORCE_X1),0)
RUST_HUF_FEATURE := huf-force-decompress-x1 RUST_HUF_FEATURE := huf-force-decompress-x1
RUST_BUILD_CONFIG := c$(ZSTD_LIB_COMPRESSION)-d$(ZSTD_LIB_DECOMPRESSION)-huf-force-decompress-x1 RUST_BUILD_CONFIG := $(RUST_MODULE_CONFIG)-huf-force-decompress-x1
endif endif
ifneq ($(RUST_HUF_FORCE_X2),0) ifneq ($(RUST_HUF_FORCE_X2),0)
RUST_HUF_FEATURE := huf-force-decompress-x2 RUST_HUF_FEATURE := huf-force-decompress-x2
RUST_BUILD_CONFIG := c$(ZSTD_LIB_COMPRESSION)-d$(ZSTD_LIB_DECOMPRESSION)-huf-force-decompress-x2 RUST_BUILD_CONFIG := $(RUST_MODULE_CONFIG)-huf-force-decompress-x2
endif endif
ifneq ($(RUST_HUF_FEATURE),) ifneq ($(RUST_HUF_FEATURE),)
ifneq ($(ZSTD_LIB_DECOMPRESSION),0) ifneq ($(ZSTD_LIB_DECOMPRESSION),0)
+5 -1886
View File
@@ -24,1890 +24,9 @@
* OTHER DEALINGS IN THE SOFTWARE. * OTHER DEALINGS IN THE SOFTWARE.
*/ */
/*- Compiler specifics -*/ /* divsufsort() is implemented in rust/src/divsufsort.rs, which provides the
#ifdef __clang__ * symbol directly. This translation unit keeps the header's prototypes in
#pragma clang diagnostic ignored "-Wshorten-64-to-32" * the build so the dictionary builder continues to compile against the
#endif * original interface. divbwt() has no callers in zstd and is declaration-
* only; it moves to Rust if a user ever appears. */
#if defined(_MSC_VER)
# pragma warning(disable : 4244)
# pragma warning(disable : 4127) /* C4127 : Condition expression is constant */
#endif
/*- Dependencies -*/
#include <assert.h>
#include <stdio.h>
#include <stdlib.h>
#include "divsufsort.h" #include "divsufsort.h"
/*- Constants -*/
#if defined(INLINE)
# undef INLINE
#endif
#if !defined(INLINE)
# define INLINE __inline
#endif
#if defined(ALPHABET_SIZE) && (ALPHABET_SIZE < 1)
# undef ALPHABET_SIZE
#endif
#if !defined(ALPHABET_SIZE)
# define ALPHABET_SIZE (256)
#endif
#define BUCKET_A_SIZE (ALPHABET_SIZE)
#define BUCKET_B_SIZE (ALPHABET_SIZE * ALPHABET_SIZE)
#if defined(SS_INSERTIONSORT_THRESHOLD)
# if SS_INSERTIONSORT_THRESHOLD < 1
# undef SS_INSERTIONSORT_THRESHOLD
# define SS_INSERTIONSORT_THRESHOLD (1)
# endif
#else
# define SS_INSERTIONSORT_THRESHOLD (8)
#endif
#if defined(SS_BLOCKSIZE)
# if SS_BLOCKSIZE < 0
# undef SS_BLOCKSIZE
# define SS_BLOCKSIZE (0)
# elif 32768 <= SS_BLOCKSIZE
# undef SS_BLOCKSIZE
# define SS_BLOCKSIZE (32767)
# endif
#else
# define SS_BLOCKSIZE (1024)
#endif
/* minstacksize = log(SS_BLOCKSIZE) / log(3) * 2 */
#if SS_BLOCKSIZE == 0
# define SS_MISORT_STACKSIZE (96)
#elif SS_BLOCKSIZE <= 4096
# define SS_MISORT_STACKSIZE (16)
#else
# define SS_MISORT_STACKSIZE (24)
#endif
#define SS_SMERGE_STACKSIZE (32)
#define TR_INSERTIONSORT_THRESHOLD (8)
#define TR_STACKSIZE (64)
/*- Macros -*/
#ifndef SWAP
# define SWAP(_a, _b) do { t = (_a); (_a) = (_b); (_b) = t; } while(0)
#endif /* SWAP */
#ifndef MIN
# define MIN(_a, _b) (((_a) < (_b)) ? (_a) : (_b))
#endif /* MIN */
#ifndef MAX
# define MAX(_a, _b) (((_a) > (_b)) ? (_a) : (_b))
#endif /* MAX */
#define STACK_PUSH(_a, _b, _c, _d)\
do {\
assert(ssize < STACK_SIZE);\
stack[ssize].a = (_a), stack[ssize].b = (_b),\
stack[ssize].c = (_c), stack[ssize++].d = (_d);\
} while(0)
#define STACK_PUSH5(_a, _b, _c, _d, _e)\
do {\
assert(ssize < STACK_SIZE);\
stack[ssize].a = (_a), stack[ssize].b = (_b),\
stack[ssize].c = (_c), stack[ssize].d = (_d), stack[ssize++].e = (_e);\
} while(0)
#define STACK_POP(_a, _b, _c, _d)\
do {\
assert(0 <= ssize);\
if(ssize == 0) { return; }\
(_a) = stack[--ssize].a, (_b) = stack[ssize].b,\
(_c) = stack[ssize].c, (_d) = stack[ssize].d;\
} while(0)
#define STACK_POP5(_a, _b, _c, _d, _e)\
do {\
assert(0 <= ssize);\
if(ssize == 0) { return; }\
(_a) = stack[--ssize].a, (_b) = stack[ssize].b,\
(_c) = stack[ssize].c, (_d) = stack[ssize].d, (_e) = stack[ssize].e;\
} while(0)
#define BUCKET_A(_c0) bucket_A[(_c0)]
#if ALPHABET_SIZE == 256
#define BUCKET_B(_c0, _c1) (bucket_B[((_c1) << 8) | (_c0)])
#define BUCKET_BSTAR(_c0, _c1) (bucket_B[((_c0) << 8) | (_c1)])
#else
#define BUCKET_B(_c0, _c1) (bucket_B[(_c1) * ALPHABET_SIZE + (_c0)])
#define BUCKET_BSTAR(_c0, _c1) (bucket_B[(_c0) * ALPHABET_SIZE + (_c1)])
#endif
/*- Private Functions -*/
static const int lg_table[256]= {
-1,0,1,1,2,2,2,2,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,
5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7
};
#if (SS_BLOCKSIZE == 0) || (SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE)
static INLINE
int
ss_ilg(int n) {
#if SS_BLOCKSIZE == 0
return (n & 0xffff0000) ?
((n & 0xff000000) ?
24 + lg_table[(n >> 24) & 0xff] :
16 + lg_table[(n >> 16) & 0xff]) :
((n & 0x0000ff00) ?
8 + lg_table[(n >> 8) & 0xff] :
0 + lg_table[(n >> 0) & 0xff]);
#elif SS_BLOCKSIZE < 256
return lg_table[n];
#else
return (n & 0xff00) ?
8 + lg_table[(n >> 8) & 0xff] :
0 + lg_table[(n >> 0) & 0xff];
#endif
}
#endif /* (SS_BLOCKSIZE == 0) || (SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE) */
#if SS_BLOCKSIZE != 0
static const int sqq_table[256] = {
0, 16, 22, 27, 32, 35, 39, 42, 45, 48, 50, 53, 55, 57, 59, 61,
64, 65, 67, 69, 71, 73, 75, 76, 78, 80, 81, 83, 84, 86, 87, 89,
90, 91, 93, 94, 96, 97, 98, 99, 101, 102, 103, 104, 106, 107, 108, 109,
110, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126,
128, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142,
143, 144, 144, 145, 146, 147, 148, 149, 150, 150, 151, 152, 153, 154, 155, 155,
156, 157, 158, 159, 160, 160, 161, 162, 163, 163, 164, 165, 166, 167, 167, 168,
169, 170, 170, 171, 172, 173, 173, 174, 175, 176, 176, 177, 178, 178, 179, 180,
181, 181, 182, 183, 183, 184, 185, 185, 186, 187, 187, 188, 189, 189, 190, 191,
192, 192, 193, 193, 194, 195, 195, 196, 197, 197, 198, 199, 199, 200, 201, 201,
202, 203, 203, 204, 204, 205, 206, 206, 207, 208, 208, 209, 209, 210, 211, 211,
212, 212, 213, 214, 214, 215, 215, 216, 217, 217, 218, 218, 219, 219, 220, 221,
221, 222, 222, 223, 224, 224, 225, 225, 226, 226, 227, 227, 228, 229, 229, 230,
230, 231, 231, 232, 232, 233, 234, 234, 235, 235, 236, 236, 237, 237, 238, 238,
239, 240, 240, 241, 241, 242, 242, 243, 243, 244, 244, 245, 245, 246, 246, 247,
247, 248, 248, 249, 249, 250, 250, 251, 251, 252, 252, 253, 253, 254, 254, 255
};
static INLINE
int
ss_isqrt(int x) {
int y, e;
if(x >= (SS_BLOCKSIZE * SS_BLOCKSIZE)) { return SS_BLOCKSIZE; }
e = (x & 0xffff0000) ?
((x & 0xff000000) ?
24 + lg_table[(x >> 24) & 0xff] :
16 + lg_table[(x >> 16) & 0xff]) :
((x & 0x0000ff00) ?
8 + lg_table[(x >> 8) & 0xff] :
0 + lg_table[(x >> 0) & 0xff]);
if(e >= 16) {
y = sqq_table[x >> ((e - 6) - (e & 1))] << ((e >> 1) - 7);
if(e >= 24) { y = (y + 1 + x / y) >> 1; }
y = (y + 1 + x / y) >> 1;
} else if(e >= 8) {
y = (sqq_table[x >> ((e - 6) - (e & 1))] >> (7 - (e >> 1))) + 1;
} else {
return sqq_table[x] >> 4;
}
return (x < (y * y)) ? y - 1 : y;
}
#endif /* SS_BLOCKSIZE != 0 */
/*---------------------------------------------------------------------------*/
/* Compares two suffixes. */
static INLINE
int
ss_compare(const unsigned char *T,
const int *p1, const int *p2,
int depth) {
const unsigned char *U1, *U2, *U1n, *U2n;
for(U1 = T + depth + *p1,
U2 = T + depth + *p2,
U1n = T + *(p1 + 1) + 2,
U2n = T + *(p2 + 1) + 2;
(U1 < U1n) && (U2 < U2n) && (*U1 == *U2);
++U1, ++U2) {
}
return U1 < U1n ?
(U2 < U2n ? *U1 - *U2 : 1) :
(U2 < U2n ? -1 : 0);
}
/*---------------------------------------------------------------------------*/
#if (SS_BLOCKSIZE != 1) && (SS_INSERTIONSORT_THRESHOLD != 1)
/* Insertionsort for small size groups */
static
void
ss_insertionsort(const unsigned char *T, const int *PA,
int *first, int *last, int depth) {
int *i, *j;
int t;
int r;
for(i = last - 2; first <= i; --i) {
for(t = *i, j = i + 1; 0 < (r = ss_compare(T, PA + t, PA + *j, depth));) {
do { *(j - 1) = *j; } while((++j < last) && (*j < 0));
if(last <= j) { break; }
}
if(r == 0) { *j = ~*j; }
*(j - 1) = t;
}
}
#endif /* (SS_BLOCKSIZE != 1) && (SS_INSERTIONSORT_THRESHOLD != 1) */
/*---------------------------------------------------------------------------*/
#if (SS_BLOCKSIZE == 0) || (SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE)
static INLINE
void
ss_fixdown(const unsigned char *Td, const int *PA,
int *SA, int i, int size) {
int j, k;
int v;
int c, d, e;
for(v = SA[i], c = Td[PA[v]]; (j = 2 * i + 1) < size; SA[i] = SA[k], i = k) {
d = Td[PA[SA[k = j++]]];
if(d < (e = Td[PA[SA[j]]])) { k = j; d = e; }
if(d <= c) { break; }
}
SA[i] = v;
}
/* Simple top-down heapsort. */
static
void
ss_heapsort(const unsigned char *Td, const int *PA, int *SA, int size) {
int i, m;
int t;
m = size;
if((size % 2) == 0) {
m--;
if(Td[PA[SA[m / 2]]] < Td[PA[SA[m]]]) { SWAP(SA[m], SA[m / 2]); }
}
for(i = m / 2 - 1; 0 <= i; --i) { ss_fixdown(Td, PA, SA, i, m); }
if((size % 2) == 0) { SWAP(SA[0], SA[m]); ss_fixdown(Td, PA, SA, 0, m); }
for(i = m - 1; 0 < i; --i) {
t = SA[0], SA[0] = SA[i];
ss_fixdown(Td, PA, SA, 0, i);
SA[i] = t;
}
}
/*---------------------------------------------------------------------------*/
/* Returns the median of three elements. */
static INLINE
int *
ss_median3(const unsigned char *Td, const int *PA,
int *v1, int *v2, int *v3) {
int *t;
if(Td[PA[*v1]] > Td[PA[*v2]]) { SWAP(v1, v2); }
if(Td[PA[*v2]] > Td[PA[*v3]]) {
if(Td[PA[*v1]] > Td[PA[*v3]]) { return v1; }
else { return v3; }
}
return v2;
}
/* Returns the median of five elements. */
static INLINE
int *
ss_median5(const unsigned char *Td, const int *PA,
int *v1, int *v2, int *v3, int *v4, int *v5) {
int *t;
if(Td[PA[*v2]] > Td[PA[*v3]]) { SWAP(v2, v3); }
if(Td[PA[*v4]] > Td[PA[*v5]]) { SWAP(v4, v5); }
if(Td[PA[*v2]] > Td[PA[*v4]]) { SWAP(v2, v4); SWAP(v3, v5); }
if(Td[PA[*v1]] > Td[PA[*v3]]) { SWAP(v1, v3); }
if(Td[PA[*v1]] > Td[PA[*v4]]) { SWAP(v1, v4); SWAP(v3, v5); }
if(Td[PA[*v3]] > Td[PA[*v4]]) { return v4; }
return v3;
}
/* Returns the pivot element. */
static INLINE
int *
ss_pivot(const unsigned char *Td, const int *PA, int *first, int *last) {
int *middle;
int t;
t = last - first;
middle = first + t / 2;
if(t <= 512) {
if(t <= 32) {
return ss_median3(Td, PA, first, middle, last - 1);
} else {
t >>= 2;
return ss_median5(Td, PA, first, first + t, middle, last - 1 - t, last - 1);
}
}
t >>= 3;
first = ss_median3(Td, PA, first, first + t, first + (t << 1));
middle = ss_median3(Td, PA, middle - t, middle, middle + t);
last = ss_median3(Td, PA, last - 1 - (t << 1), last - 1 - t, last - 1);
return ss_median3(Td, PA, first, middle, last);
}
/*---------------------------------------------------------------------------*/
/* Binary partition for substrings. */
static INLINE
int *
ss_partition(const int *PA,
int *first, int *last, int depth) {
int *a, *b;
int t;
for(a = first - 1, b = last;;) {
for(; (++a < b) && ((PA[*a] + depth) >= (PA[*a + 1] + 1));) { *a = ~*a; }
for(; (a < --b) && ((PA[*b] + depth) < (PA[*b + 1] + 1));) { }
if(b <= a) { break; }
t = ~*b;
*b = *a;
*a = t;
}
if(first < a) { *first = ~*first; }
return a;
}
/* Multikey introsort for medium size groups. */
static
void
ss_mintrosort(const unsigned char *T, const int *PA,
int *first, int *last,
int depth) {
#define STACK_SIZE SS_MISORT_STACKSIZE
struct { int *a, *b, c; int d; } stack[STACK_SIZE];
const unsigned char *Td;
int *a, *b, *c, *d, *e, *f;
int s, t;
int ssize;
int limit;
int v, x = 0;
for(ssize = 0, limit = ss_ilg(last - first);;) {
if((last - first) <= SS_INSERTIONSORT_THRESHOLD) {
#if 1 < SS_INSERTIONSORT_THRESHOLD
if(1 < (last - first)) { ss_insertionsort(T, PA, first, last, depth); }
#endif
STACK_POP(first, last, depth, limit);
continue;
}
Td = T + depth;
if(limit-- == 0) { ss_heapsort(Td, PA, first, last - first); }
if(limit < 0) {
for(a = first + 1, v = Td[PA[*first]]; a < last; ++a) {
if((x = Td[PA[*a]]) != v) {
if(1 < (a - first)) { break; }
v = x;
first = a;
}
}
if(Td[PA[*first] - 1] < v) {
first = ss_partition(PA, first, a, depth);
}
if((a - first) <= (last - a)) {
if(1 < (a - first)) {
STACK_PUSH(a, last, depth, -1);
last = a, depth += 1, limit = ss_ilg(a - first);
} else {
first = a, limit = -1;
}
} else {
if(1 < (last - a)) {
STACK_PUSH(first, a, depth + 1, ss_ilg(a - first));
first = a, limit = -1;
} else {
last = a, depth += 1, limit = ss_ilg(a - first);
}
}
continue;
}
/* choose pivot */
a = ss_pivot(Td, PA, first, last);
v = Td[PA[*a]];
SWAP(*first, *a);
/* partition */
for(b = first; (++b < last) && ((x = Td[PA[*b]]) == v);) { }
if(((a = b) < last) && (x < v)) {
for(; (++b < last) && ((x = Td[PA[*b]]) <= v);) {
if(x == v) { SWAP(*b, *a); ++a; }
}
}
for(c = last; (b < --c) && ((x = Td[PA[*c]]) == v);) { }
if((b < (d = c)) && (x > v)) {
for(; (b < --c) && ((x = Td[PA[*c]]) >= v);) {
if(x == v) { SWAP(*c, *d); --d; }
}
}
for(; b < c;) {
SWAP(*b, *c);
for(; (++b < c) && ((x = Td[PA[*b]]) <= v);) {
if(x == v) { SWAP(*b, *a); ++a; }
}
for(; (b < --c) && ((x = Td[PA[*c]]) >= v);) {
if(x == v) { SWAP(*c, *d); --d; }
}
}
if(a <= d) {
c = b - 1;
if((s = a - first) > (t = b - a)) { s = t; }
for(e = first, f = b - s; 0 < s; --s, ++e, ++f) { SWAP(*e, *f); }
if((s = d - c) > (t = last - d - 1)) { s = t; }
for(e = b, f = last - s; 0 < s; --s, ++e, ++f) { SWAP(*e, *f); }
a = first + (b - a), c = last - (d - c);
b = (v <= Td[PA[*a] - 1]) ? a : ss_partition(PA, a, c, depth);
if((a - first) <= (last - c)) {
if((last - c) <= (c - b)) {
STACK_PUSH(b, c, depth + 1, ss_ilg(c - b));
STACK_PUSH(c, last, depth, limit);
last = a;
} else if((a - first) <= (c - b)) {
STACK_PUSH(c, last, depth, limit);
STACK_PUSH(b, c, depth + 1, ss_ilg(c - b));
last = a;
} else {
STACK_PUSH(c, last, depth, limit);
STACK_PUSH(first, a, depth, limit);
first = b, last = c, depth += 1, limit = ss_ilg(c - b);
}
} else {
if((a - first) <= (c - b)) {
STACK_PUSH(b, c, depth + 1, ss_ilg(c - b));
STACK_PUSH(first, a, depth, limit);
first = c;
} else if((last - c) <= (c - b)) {
STACK_PUSH(first, a, depth, limit);
STACK_PUSH(b, c, depth + 1, ss_ilg(c - b));
first = c;
} else {
STACK_PUSH(first, a, depth, limit);
STACK_PUSH(c, last, depth, limit);
first = b, last = c, depth += 1, limit = ss_ilg(c - b);
}
}
} else {
limit += 1;
if(Td[PA[*first] - 1] < v) {
first = ss_partition(PA, first, last, depth);
limit = ss_ilg(last - first);
}
depth += 1;
}
}
#undef STACK_SIZE
}
#endif /* (SS_BLOCKSIZE == 0) || (SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE) */
/*---------------------------------------------------------------------------*/
#if SS_BLOCKSIZE != 0
static INLINE
void
ss_blockswap(int *a, int *b, int n) {
int t;
for(; 0 < n; --n, ++a, ++b) {
t = *a, *a = *b, *b = t;
}
}
static INLINE
void
ss_rotate(int *first, int *middle, int *last) {
int *a, *b, t;
int l, r;
l = middle - first, r = last - middle;
for(; (0 < l) && (0 < r);) {
if(l == r) { ss_blockswap(first, middle, l); break; }
if(l < r) {
a = last - 1, b = middle - 1;
t = *a;
do {
*a-- = *b, *b-- = *a;
if(b < first) {
*a = t;
last = a;
if((r -= l + 1) <= l) { break; }
a -= 1, b = middle - 1;
t = *a;
}
} while(1);
} else {
a = first, b = middle;
t = *a;
do {
*a++ = *b, *b++ = *a;
if(last <= b) {
*a = t;
first = a + 1;
if((l -= r + 1) <= r) { break; }
a += 1, b = middle;
t = *a;
}
} while(1);
}
}
}
/*---------------------------------------------------------------------------*/
static
void
ss_inplacemerge(const unsigned char *T, const int *PA,
int *first, int *middle, int *last,
int depth) {
const int *p;
int *a, *b;
int len, half;
int q, r;
int x;
for(;;) {
if(*(last - 1) < 0) { x = 1; p = PA + ~*(last - 1); }
else { x = 0; p = PA + *(last - 1); }
for(a = first, len = middle - first, half = len >> 1, r = -1;
0 < len;
len = half, half >>= 1) {
b = a + half;
q = ss_compare(T, PA + ((0 <= *b) ? *b : ~*b), p, depth);
if(q < 0) {
a = b + 1;
half -= (len & 1) ^ 1;
} else {
r = q;
}
}
if(a < middle) {
if(r == 0) { *a = ~*a; }
ss_rotate(a, middle, last);
last -= middle - a;
middle = a;
if(first == middle) { break; }
}
--last;
if(x != 0) { while(*--last < 0) { } }
if(middle == last) { break; }
}
}
/*---------------------------------------------------------------------------*/
/* Merge-forward with internal buffer. */
static
void
ss_mergeforward(const unsigned char *T, const int *PA,
int *first, int *middle, int *last,
int *buf, int depth) {
int *a, *b, *c, *bufend;
int t;
int r;
bufend = buf + (middle - first) - 1;
ss_blockswap(buf, first, middle - first);
for(t = *(a = first), b = buf, c = middle;;) {
r = ss_compare(T, PA + *b, PA + *c, depth);
if(r < 0) {
do {
*a++ = *b;
if(bufend <= b) { *bufend = t; return; }
*b++ = *a;
} while(*b < 0);
} else if(r > 0) {
do {
*a++ = *c, *c++ = *a;
if(last <= c) {
while(b < bufend) { *a++ = *b, *b++ = *a; }
*a = *b, *b = t;
return;
}
} while(*c < 0);
} else {
*c = ~*c;
do {
*a++ = *b;
if(bufend <= b) { *bufend = t; return; }
*b++ = *a;
} while(*b < 0);
do {
*a++ = *c, *c++ = *a;
if(last <= c) {
while(b < bufend) { *a++ = *b, *b++ = *a; }
*a = *b, *b = t;
return;
}
} while(*c < 0);
}
}
}
/* Merge-backward with internal buffer. */
static
void
ss_mergebackward(const unsigned char *T, const int *PA,
int *first, int *middle, int *last,
int *buf, int depth) {
const int *p1, *p2;
int *a, *b, *c, *bufend;
int t;
int r;
int x;
bufend = buf + (last - middle) - 1;
ss_blockswap(buf, middle, last - middle);
x = 0;
if(*bufend < 0) { p1 = PA + ~*bufend; x |= 1; }
else { p1 = PA + *bufend; }
if(*(middle - 1) < 0) { p2 = PA + ~*(middle - 1); x |= 2; }
else { p2 = PA + *(middle - 1); }
for(t = *(a = last - 1), b = bufend, c = middle - 1;;) {
r = ss_compare(T, p1, p2, depth);
if(0 < r) {
if(x & 1) { do { *a-- = *b, *b-- = *a; } while(*b < 0); x ^= 1; }
*a-- = *b;
if(b <= buf) { *buf = t; break; }
*b-- = *a;
if(*b < 0) { p1 = PA + ~*b; x |= 1; }
else { p1 = PA + *b; }
} else if(r < 0) {
if(x & 2) { do { *a-- = *c, *c-- = *a; } while(*c < 0); x ^= 2; }
*a-- = *c, *c-- = *a;
if(c < first) {
while(buf < b) { *a-- = *b, *b-- = *a; }
*a = *b, *b = t;
break;
}
if(*c < 0) { p2 = PA + ~*c; x |= 2; }
else { p2 = PA + *c; }
} else {
if(x & 1) { do { *a-- = *b, *b-- = *a; } while(*b < 0); x ^= 1; }
*a-- = ~*b;
if(b <= buf) { *buf = t; break; }
*b-- = *a;
if(x & 2) { do { *a-- = *c, *c-- = *a; } while(*c < 0); x ^= 2; }
*a-- = *c, *c-- = *a;
if(c < first) {
while(buf < b) { *a-- = *b, *b-- = *a; }
*a = *b, *b = t;
break;
}
if(*b < 0) { p1 = PA + ~*b; x |= 1; }
else { p1 = PA + *b; }
if(*c < 0) { p2 = PA + ~*c; x |= 2; }
else { p2 = PA + *c; }
}
}
}
/* D&C based merge. */
static
void
ss_swapmerge(const unsigned char *T, const int *PA,
int *first, int *middle, int *last,
int *buf, int bufsize, int depth) {
#define STACK_SIZE SS_SMERGE_STACKSIZE
#define GETIDX(a) ((0 <= (a)) ? (a) : (~(a)))
#define MERGE_CHECK(a, b, c)\
do {\
if(((c) & 1) ||\
(((c) & 2) && (ss_compare(T, PA + GETIDX(*((a) - 1)), PA + *(a), depth) == 0))) {\
*(a) = ~*(a);\
}\
if(((c) & 4) && ((ss_compare(T, PA + GETIDX(*((b) - 1)), PA + *(b), depth) == 0))) {\
*(b) = ~*(b);\
}\
} while(0)
struct { int *a, *b, *c; int d; } stack[STACK_SIZE];
int *l, *r, *lm, *rm;
int m, len, half;
int ssize;
int check, next;
for(check = 0, ssize = 0;;) {
if((last - middle) <= bufsize) {
if((first < middle) && (middle < last)) {
ss_mergebackward(T, PA, first, middle, last, buf, depth);
}
MERGE_CHECK(first, last, check);
STACK_POP(first, middle, last, check);
continue;
}
if((middle - first) <= bufsize) {
if(first < middle) {
ss_mergeforward(T, PA, first, middle, last, buf, depth);
}
MERGE_CHECK(first, last, check);
STACK_POP(first, middle, last, check);
continue;
}
for(m = 0, len = MIN(middle - first, last - middle), half = len >> 1;
0 < len;
len = half, half >>= 1) {
if(ss_compare(T, PA + GETIDX(*(middle + m + half)),
PA + GETIDX(*(middle - m - half - 1)), depth) < 0) {
m += half + 1;
half -= (len & 1) ^ 1;
}
}
if(0 < m) {
lm = middle - m, rm = middle + m;
ss_blockswap(lm, middle, m);
l = r = middle, next = 0;
if(rm < last) {
if(*rm < 0) {
*rm = ~*rm;
if(first < lm) { for(; *--l < 0;) { } next |= 4; }
next |= 1;
} else if(first < lm) {
for(; *r < 0; ++r) { }
next |= 2;
}
}
if((l - first) <= (last - r)) {
STACK_PUSH(r, rm, last, (next & 3) | (check & 4));
middle = lm, last = l, check = (check & 3) | (next & 4);
} else {
if((next & 2) && (r == middle)) { next ^= 6; }
STACK_PUSH(first, lm, l, (check & 3) | (next & 4));
first = r, middle = rm, check = (next & 3) | (check & 4);
}
} else {
if(ss_compare(T, PA + GETIDX(*(middle - 1)), PA + *middle, depth) == 0) {
*middle = ~*middle;
}
MERGE_CHECK(first, last, check);
STACK_POP(first, middle, last, check);
}
}
#undef STACK_SIZE
}
#endif /* SS_BLOCKSIZE != 0 */
/*---------------------------------------------------------------------------*/
/* Substring sort */
static
void
sssort(const unsigned char *T, const int *PA,
int *first, int *last,
int *buf, int bufsize,
int depth, int n, int lastsuffix) {
int *a;
#if SS_BLOCKSIZE != 0
int *b, *middle, *curbuf;
int j, k, curbufsize, limit;
#endif
int i;
if(lastsuffix != 0) { ++first; }
#if SS_BLOCKSIZE == 0
ss_mintrosort(T, PA, first, last, depth);
#else
if((bufsize < SS_BLOCKSIZE) &&
(bufsize < (last - first)) &&
(bufsize < (limit = ss_isqrt(last - first)))) {
if(SS_BLOCKSIZE < limit) { limit = SS_BLOCKSIZE; }
buf = middle = last - limit, bufsize = limit;
} else {
middle = last, limit = 0;
}
for(a = first, i = 0; SS_BLOCKSIZE < (middle - a); a += SS_BLOCKSIZE, ++i) {
#if SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE
ss_mintrosort(T, PA, a, a + SS_BLOCKSIZE, depth);
#elif 1 < SS_BLOCKSIZE
ss_insertionsort(T, PA, a, a + SS_BLOCKSIZE, depth);
#endif
curbufsize = last - (a + SS_BLOCKSIZE);
curbuf = a + SS_BLOCKSIZE;
if(curbufsize <= bufsize) { curbufsize = bufsize, curbuf = buf; }
for(b = a, k = SS_BLOCKSIZE, j = i; j & 1; b -= k, k <<= 1, j >>= 1) {
ss_swapmerge(T, PA, b - k, b, b + k, curbuf, curbufsize, depth);
}
}
#if SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE
ss_mintrosort(T, PA, a, middle, depth);
#elif 1 < SS_BLOCKSIZE
ss_insertionsort(T, PA, a, middle, depth);
#endif
for(k = SS_BLOCKSIZE; i != 0; k <<= 1, i >>= 1) {
if(i & 1) {
ss_swapmerge(T, PA, a - k, a, middle, buf, bufsize, depth);
a -= k;
}
}
if(limit != 0) {
#if SS_INSERTIONSORT_THRESHOLD < SS_BLOCKSIZE
ss_mintrosort(T, PA, middle, last, depth);
#elif 1 < SS_BLOCKSIZE
ss_insertionsort(T, PA, middle, last, depth);
#endif
ss_inplacemerge(T, PA, first, middle, last, depth);
}
#endif
if(lastsuffix != 0) {
/* Insert last type B* suffix. */
int PAi[2]; PAi[0] = PA[*(first - 1)], PAi[1] = n - 2;
for(a = first, i = *(first - 1);
(a < last) && ((*a < 0) || (0 < ss_compare(T, &(PAi[0]), PA + *a, depth)));
++a) {
*(a - 1) = *a;
}
*(a - 1) = i;
}
}
/*---------------------------------------------------------------------------*/
static INLINE
int
tr_ilg(int n) {
return (n & 0xffff0000) ?
((n & 0xff000000) ?
24 + lg_table[(n >> 24) & 0xff] :
16 + lg_table[(n >> 16) & 0xff]) :
((n & 0x0000ff00) ?
8 + lg_table[(n >> 8) & 0xff] :
0 + lg_table[(n >> 0) & 0xff]);
}
/*---------------------------------------------------------------------------*/
/* Simple insertionsort for small size groups. */
static
void
tr_insertionsort(const int *ISAd, int *first, int *last) {
int *a, *b;
int t, r;
for(a = first + 1; a < last; ++a) {
for(t = *a, b = a - 1; 0 > (r = ISAd[t] - ISAd[*b]);) {
do { *(b + 1) = *b; } while((first <= --b) && (*b < 0));
if(b < first) { break; }
}
if(r == 0) { *b = ~*b; }
*(b + 1) = t;
}
}
/*---------------------------------------------------------------------------*/
static INLINE
void
tr_fixdown(const int *ISAd, int *SA, int i, int size) {
int j, k;
int v;
int c, d, e;
for(v = SA[i], c = ISAd[v]; (j = 2 * i + 1) < size; SA[i] = SA[k], i = k) {
d = ISAd[SA[k = j++]];
if(d < (e = ISAd[SA[j]])) { k = j; d = e; }
if(d <= c) { break; }
}
SA[i] = v;
}
/* Simple top-down heapsort. */
static
void
tr_heapsort(const int *ISAd, int *SA, int size) {
int i, m;
int t;
m = size;
if((size % 2) == 0) {
m--;
if(ISAd[SA[m / 2]] < ISAd[SA[m]]) { SWAP(SA[m], SA[m / 2]); }
}
for(i = m / 2 - 1; 0 <= i; --i) { tr_fixdown(ISAd, SA, i, m); }
if((size % 2) == 0) { SWAP(SA[0], SA[m]); tr_fixdown(ISAd, SA, 0, m); }
for(i = m - 1; 0 < i; --i) {
t = SA[0], SA[0] = SA[i];
tr_fixdown(ISAd, SA, 0, i);
SA[i] = t;
}
}
/*---------------------------------------------------------------------------*/
/* Returns the median of three elements. */
static INLINE
int *
tr_median3(const int *ISAd, int *v1, int *v2, int *v3) {
int *t;
if(ISAd[*v1] > ISAd[*v2]) { SWAP(v1, v2); }
if(ISAd[*v2] > ISAd[*v3]) {
if(ISAd[*v1] > ISAd[*v3]) { return v1; }
else { return v3; }
}
return v2;
}
/* Returns the median of five elements. */
static INLINE
int *
tr_median5(const int *ISAd,
int *v1, int *v2, int *v3, int *v4, int *v5) {
int *t;
if(ISAd[*v2] > ISAd[*v3]) { SWAP(v2, v3); }
if(ISAd[*v4] > ISAd[*v5]) { SWAP(v4, v5); }
if(ISAd[*v2] > ISAd[*v4]) { SWAP(v2, v4); SWAP(v3, v5); }
if(ISAd[*v1] > ISAd[*v3]) { SWAP(v1, v3); }
if(ISAd[*v1] > ISAd[*v4]) { SWAP(v1, v4); SWAP(v3, v5); }
if(ISAd[*v3] > ISAd[*v4]) { return v4; }
return v3;
}
/* Returns the pivot element. */
static INLINE
int *
tr_pivot(const int *ISAd, int *first, int *last) {
int *middle;
int t;
t = last - first;
middle = first + t / 2;
if(t <= 512) {
if(t <= 32) {
return tr_median3(ISAd, first, middle, last - 1);
} else {
t >>= 2;
return tr_median5(ISAd, first, first + t, middle, last - 1 - t, last - 1);
}
}
t >>= 3;
first = tr_median3(ISAd, first, first + t, first + (t << 1));
middle = tr_median3(ISAd, middle - t, middle, middle + t);
last = tr_median3(ISAd, last - 1 - (t << 1), last - 1 - t, last - 1);
return tr_median3(ISAd, first, middle, last);
}
/*---------------------------------------------------------------------------*/
typedef struct _trbudget_t trbudget_t;
struct _trbudget_t {
int chance;
int remain;
int incval;
int count;
};
static INLINE
void
trbudget_init(trbudget_t *budget, int chance, int incval) {
budget->chance = chance;
budget->remain = budget->incval = incval;
}
static INLINE
int
trbudget_check(trbudget_t *budget, int size) {
if(size <= budget->remain) { budget->remain -= size; return 1; }
if(budget->chance == 0) { budget->count += size; return 0; }
budget->remain += budget->incval - size;
budget->chance -= 1;
return 1;
}
/*---------------------------------------------------------------------------*/
static INLINE
void
tr_partition(const int *ISAd,
int *first, int *middle, int *last,
int **pa, int **pb, int v) {
int *a, *b, *c, *d, *e, *f;
int t, s;
int x = 0;
for(b = middle - 1; (++b < last) && ((x = ISAd[*b]) == v);) { }
if(((a = b) < last) && (x < v)) {
for(; (++b < last) && ((x = ISAd[*b]) <= v);) {
if(x == v) { SWAP(*b, *a); ++a; }
}
}
for(c = last; (b < --c) && ((x = ISAd[*c]) == v);) { }
if((b < (d = c)) && (x > v)) {
for(; (b < --c) && ((x = ISAd[*c]) >= v);) {
if(x == v) { SWAP(*c, *d); --d; }
}
}
for(; b < c;) {
SWAP(*b, *c);
for(; (++b < c) && ((x = ISAd[*b]) <= v);) {
if(x == v) { SWAP(*b, *a); ++a; }
}
for(; (b < --c) && ((x = ISAd[*c]) >= v);) {
if(x == v) { SWAP(*c, *d); --d; }
}
}
if(a <= d) {
c = b - 1;
if((s = a - first) > (t = b - a)) { s = t; }
for(e = first, f = b - s; 0 < s; --s, ++e, ++f) { SWAP(*e, *f); }
if((s = d - c) > (t = last - d - 1)) { s = t; }
for(e = b, f = last - s; 0 < s; --s, ++e, ++f) { SWAP(*e, *f); }
first += (b - a), last -= (d - c);
}
*pa = first, *pb = last;
}
static
void
tr_copy(int *ISA, const int *SA,
int *first, int *a, int *b, int *last,
int depth) {
/* sort suffixes of middle partition
by using sorted order of suffixes of left and right partition. */
int *c, *d, *e;
int s, v;
v = b - SA - 1;
for(c = first, d = a - 1; c <= d; ++c) {
if((0 <= (s = *c - depth)) && (ISA[s] == v)) {
*++d = s;
ISA[s] = d - SA;
}
}
for(c = last - 1, e = d + 1, d = b; e < d; --c) {
if((0 <= (s = *c - depth)) && (ISA[s] == v)) {
*--d = s;
ISA[s] = d - SA;
}
}
}
static
void
tr_partialcopy(int *ISA, const int *SA,
int *first, int *a, int *b, int *last,
int depth) {
int *c, *d, *e;
int s, v;
int rank, lastrank, newrank = -1;
v = b - SA - 1;
lastrank = -1;
for(c = first, d = a - 1; c <= d; ++c) {
if((0 <= (s = *c - depth)) && (ISA[s] == v)) {
*++d = s;
rank = ISA[s + depth];
if(lastrank != rank) { lastrank = rank; newrank = d - SA; }
ISA[s] = newrank;
}
}
lastrank = -1;
for(e = d; first <= e; --e) {
rank = ISA[*e];
if(lastrank != rank) { lastrank = rank; newrank = e - SA; }
if(newrank != rank) { ISA[*e] = newrank; }
}
lastrank = -1;
for(c = last - 1, e = d + 1, d = b; e < d; --c) {
if((0 <= (s = *c - depth)) && (ISA[s] == v)) {
*--d = s;
rank = ISA[s + depth];
if(lastrank != rank) { lastrank = rank; newrank = d - SA; }
ISA[s] = newrank;
}
}
}
static
void
tr_introsort(int *ISA, const int *ISAd,
int *SA, int *first, int *last,
trbudget_t *budget) {
#define STACK_SIZE TR_STACKSIZE
struct { const int *a; int *b, *c; int d, e; }stack[STACK_SIZE];
int *a, *b, *c;
int t;
int v, x = 0;
int incr = ISAd - ISA;
int limit, next;
int ssize, trlink = -1;
for(ssize = 0, limit = tr_ilg(last - first);;) {
if(limit < 0) {
if(limit == -1) {
/* tandem repeat partition */
tr_partition(ISAd - incr, first, first, last, &a, &b, last - SA - 1);
/* update ranks */
if(a < last) {
for(c = first, v = a - SA - 1; c < a; ++c) { ISA[*c] = v; }
}
if(b < last) {
for(c = a, v = b - SA - 1; c < b; ++c) { ISA[*c] = v; }
}
/* push */
if(1 < (b - a)) {
STACK_PUSH5(NULL, a, b, 0, 0);
STACK_PUSH5(ISAd - incr, first, last, -2, trlink);
trlink = ssize - 2;
}
if((a - first) <= (last - b)) {
if(1 < (a - first)) {
STACK_PUSH5(ISAd, b, last, tr_ilg(last - b), trlink);
last = a, limit = tr_ilg(a - first);
} else if(1 < (last - b)) {
first = b, limit = tr_ilg(last - b);
} else {
STACK_POP5(ISAd, first, last, limit, trlink);
}
} else {
if(1 < (last - b)) {
STACK_PUSH5(ISAd, first, a, tr_ilg(a - first), trlink);
first = b, limit = tr_ilg(last - b);
} else if(1 < (a - first)) {
last = a, limit = tr_ilg(a - first);
} else {
STACK_POP5(ISAd, first, last, limit, trlink);
}
}
} else if(limit == -2) {
/* tandem repeat copy */
a = stack[--ssize].b, b = stack[ssize].c;
if(stack[ssize].d == 0) {
tr_copy(ISA, SA, first, a, b, last, ISAd - ISA);
} else {
if(0 <= trlink) { stack[trlink].d = -1; }
tr_partialcopy(ISA, SA, first, a, b, last, ISAd - ISA);
}
STACK_POP5(ISAd, first, last, limit, trlink);
} else {
/* sorted partition */
if(0 <= *first) {
a = first;
do { ISA[*a] = a - SA; } while((++a < last) && (0 <= *a));
first = a;
}
if(first < last) {
a = first; do { *a = ~*a; } while(*++a < 0);
next = (ISA[*a] != ISAd[*a]) ? tr_ilg(a - first + 1) : -1;
if(++a < last) { for(b = first, v = a - SA - 1; b < a; ++b) { ISA[*b] = v; } }
/* push */
if(trbudget_check(budget, a - first)) {
if((a - first) <= (last - a)) {
STACK_PUSH5(ISAd, a, last, -3, trlink);
ISAd += incr, last = a, limit = next;
} else {
if(1 < (last - a)) {
STACK_PUSH5(ISAd + incr, first, a, next, trlink);
first = a, limit = -3;
} else {
ISAd += incr, last = a, limit = next;
}
}
} else {
if(0 <= trlink) { stack[trlink].d = -1; }
if(1 < (last - a)) {
first = a, limit = -3;
} else {
STACK_POP5(ISAd, first, last, limit, trlink);
}
}
} else {
STACK_POP5(ISAd, first, last, limit, trlink);
}
}
continue;
}
if((last - first) <= TR_INSERTIONSORT_THRESHOLD) {
tr_insertionsort(ISAd, first, last);
limit = -3;
continue;
}
if(limit-- == 0) {
tr_heapsort(ISAd, first, last - first);
for(a = last - 1; first < a; a = b) {
for(x = ISAd[*a], b = a - 1; (first <= b) && (ISAd[*b] == x); --b) { *b = ~*b; }
}
limit = -3;
continue;
}
/* choose pivot */
a = tr_pivot(ISAd, first, last);
SWAP(*first, *a);
v = ISAd[*first];
/* partition */
tr_partition(ISAd, first, first + 1, last, &a, &b, v);
if((last - first) != (b - a)) {
next = (ISA[*a] != v) ? tr_ilg(b - a) : -1;
/* update ranks */
for(c = first, v = a - SA - 1; c < a; ++c) { ISA[*c] = v; }
if(b < last) { for(c = a, v = b - SA - 1; c < b; ++c) { ISA[*c] = v; } }
/* push */
if((1 < (b - a)) && (trbudget_check(budget, b - a))) {
if((a - first) <= (last - b)) {
if((last - b) <= (b - a)) {
if(1 < (a - first)) {
STACK_PUSH5(ISAd + incr, a, b, next, trlink);
STACK_PUSH5(ISAd, b, last, limit, trlink);
last = a;
} else if(1 < (last - b)) {
STACK_PUSH5(ISAd + incr, a, b, next, trlink);
first = b;
} else {
ISAd += incr, first = a, last = b, limit = next;
}
} else if((a - first) <= (b - a)) {
if(1 < (a - first)) {
STACK_PUSH5(ISAd, b, last, limit, trlink);
STACK_PUSH5(ISAd + incr, a, b, next, trlink);
last = a;
} else {
STACK_PUSH5(ISAd, b, last, limit, trlink);
ISAd += incr, first = a, last = b, limit = next;
}
} else {
STACK_PUSH5(ISAd, b, last, limit, trlink);
STACK_PUSH5(ISAd, first, a, limit, trlink);
ISAd += incr, first = a, last = b, limit = next;
}
} else {
if((a - first) <= (b - a)) {
if(1 < (last - b)) {
STACK_PUSH5(ISAd + incr, a, b, next, trlink);
STACK_PUSH5(ISAd, first, a, limit, trlink);
first = b;
} else if(1 < (a - first)) {
STACK_PUSH5(ISAd + incr, a, b, next, trlink);
last = a;
} else {
ISAd += incr, first = a, last = b, limit = next;
}
} else if((last - b) <= (b - a)) {
if(1 < (last - b)) {
STACK_PUSH5(ISAd, first, a, limit, trlink);
STACK_PUSH5(ISAd + incr, a, b, next, trlink);
first = b;
} else {
STACK_PUSH5(ISAd, first, a, limit, trlink);
ISAd += incr, first = a, last = b, limit = next;
}
} else {
STACK_PUSH5(ISAd, first, a, limit, trlink);
STACK_PUSH5(ISAd, b, last, limit, trlink);
ISAd += incr, first = a, last = b, limit = next;
}
}
} else {
if((1 < (b - a)) && (0 <= trlink)) { stack[trlink].d = -1; }
if((a - first) <= (last - b)) {
if(1 < (a - first)) {
STACK_PUSH5(ISAd, b, last, limit, trlink);
last = a;
} else if(1 < (last - b)) {
first = b;
} else {
STACK_POP5(ISAd, first, last, limit, trlink);
}
} else {
if(1 < (last - b)) {
STACK_PUSH5(ISAd, first, a, limit, trlink);
first = b;
} else if(1 < (a - first)) {
last = a;
} else {
STACK_POP5(ISAd, first, last, limit, trlink);
}
}
}
} else {
if(trbudget_check(budget, last - first)) {
limit = tr_ilg(last - first), ISAd += incr;
} else {
if(0 <= trlink) { stack[trlink].d = -1; }
STACK_POP5(ISAd, first, last, limit, trlink);
}
}
}
#undef STACK_SIZE
}
/*---------------------------------------------------------------------------*/
/* Tandem repeat sort */
static
void
trsort(int *ISA, int *SA, int n, int depth) {
int *ISAd;
int *first, *last;
trbudget_t budget;
int t, skip, unsorted;
trbudget_init(&budget, tr_ilg(n) * 2 / 3, n);
/* trbudget_init(&budget, tr_ilg(n) * 3 / 4, n); */
for(ISAd = ISA + depth; -n < *SA; ISAd += ISAd - ISA) {
first = SA;
skip = 0;
unsorted = 0;
do {
if((t = *first) < 0) { first -= t; skip += t; }
else {
if(skip != 0) { *(first + skip) = skip; skip = 0; }
last = SA + ISA[t] + 1;
if(1 < (last - first)) {
budget.count = 0;
tr_introsort(ISA, ISAd, SA, first, last, &budget);
if(budget.count != 0) { unsorted += budget.count; }
else { skip = first - last; }
} else if((last - first) == 1) {
skip = -1;
}
first = last;
}
} while(first < (SA + n));
if(skip != 0) { *(first + skip) = skip; }
if(unsorted == 0) { break; }
}
}
/*---------------------------------------------------------------------------*/
/* Sorts suffixes of type B*. */
static
int
sort_typeBstar(const unsigned char *T, int *SA,
int *bucket_A, int *bucket_B,
int n, int openMP) {
int *PAb, *ISAb, *buf;
#ifdef LIBBSC_OPENMP
int *curbuf;
int l;
#endif
int i, j, k, t, m, bufsize;
int c0, c1;
#ifdef LIBBSC_OPENMP
int d0, d1;
#endif
(void)openMP;
/* Initialize bucket arrays. */
for(i = 0; i < BUCKET_A_SIZE; ++i) { bucket_A[i] = 0; }
for(i = 0; i < BUCKET_B_SIZE; ++i) { bucket_B[i] = 0; }
/* Count the number of occurrences of the first one or two characters of each
type A, B and B* suffix. Moreover, store the beginning position of all
type B* suffixes into the array SA. */
for(i = n - 1, m = n, c0 = T[n - 1]; 0 <= i;) {
/* type A suffix. */
do { ++BUCKET_A(c1 = c0); } while((0 <= --i) && ((c0 = T[i]) >= c1));
if(0 <= i) {
/* type B* suffix. */
++BUCKET_BSTAR(c0, c1);
SA[--m] = i;
/* type B suffix. */
for(--i, c1 = c0; (0 <= i) && ((c0 = T[i]) <= c1); --i, c1 = c0) {
++BUCKET_B(c0, c1);
}
}
}
m = n - m;
/*
note:
A type B* suffix is lexicographically smaller than a type B suffix that
begins with the same first two characters.
*/
/* Calculate the index of start/end point of each bucket. */
for(c0 = 0, i = 0, j = 0; c0 < ALPHABET_SIZE; ++c0) {
t = i + BUCKET_A(c0);
BUCKET_A(c0) = i + j; /* start point */
i = t + BUCKET_B(c0, c0);
for(c1 = c0 + 1; c1 < ALPHABET_SIZE; ++c1) {
j += BUCKET_BSTAR(c0, c1);
BUCKET_BSTAR(c0, c1) = j; /* end point */
i += BUCKET_B(c0, c1);
}
}
if(0 < m) {
/* Sort the type B* suffixes by their first two characters. */
PAb = SA + n - m; ISAb = SA + m;
for(i = m - 2; 0 <= i; --i) {
t = PAb[i], c0 = T[t], c1 = T[t + 1];
SA[--BUCKET_BSTAR(c0, c1)] = i;
}
t = PAb[m - 1], c0 = T[t], c1 = T[t + 1];
SA[--BUCKET_BSTAR(c0, c1)] = m - 1;
/* Sort the type B* substrings using sssort. */
#ifdef LIBBSC_OPENMP
if (openMP)
{
buf = SA + m;
c0 = ALPHABET_SIZE - 2, c1 = ALPHABET_SIZE - 1, j = m;
#pragma omp parallel default(shared) private(bufsize, curbuf, k, l, d0, d1)
{
bufsize = (n - (2 * m)) / omp_get_num_threads();
curbuf = buf + omp_get_thread_num() * bufsize;
k = 0;
for(;;) {
#pragma omp critical(sssort_lock)
{
if(0 < (l = j)) {
d0 = c0, d1 = c1;
do {
k = BUCKET_BSTAR(d0, d1);
if(--d1 <= d0) {
d1 = ALPHABET_SIZE - 1;
if(--d0 < 0) { break; }
}
} while(((l - k) <= 1) && (0 < (l = k)));
c0 = d0, c1 = d1, j = k;
}
}
if(l == 0) { break; }
sssort(T, PAb, SA + k, SA + l,
curbuf, bufsize, 2, n, *(SA + k) == (m - 1));
}
}
}
else
{
buf = SA + m, bufsize = n - (2 * m);
for(c0 = ALPHABET_SIZE - 2, j = m; 0 < j; --c0) {
for(c1 = ALPHABET_SIZE - 1; c0 < c1; j = i, --c1) {
i = BUCKET_BSTAR(c0, c1);
if(1 < (j - i)) {
sssort(T, PAb, SA + i, SA + j,
buf, bufsize, 2, n, *(SA + i) == (m - 1));
}
}
}
}
#else
buf = SA + m, bufsize = n - (2 * m);
for(c0 = ALPHABET_SIZE - 2, j = m; 0 < j; --c0) {
for(c1 = ALPHABET_SIZE - 1; c0 < c1; j = i, --c1) {
i = BUCKET_BSTAR(c0, c1);
if(1 < (j - i)) {
sssort(T, PAb, SA + i, SA + j,
buf, bufsize, 2, n, *(SA + i) == (m - 1));
}
}
}
#endif
/* Compute ranks of type B* substrings. */
for(i = m - 1; 0 <= i; --i) {
if(0 <= SA[i]) {
j = i;
do { ISAb[SA[i]] = i; } while((0 <= --i) && (0 <= SA[i]));
SA[i + 1] = i - j;
if(i <= 0) { break; }
}
j = i;
do { ISAb[SA[i] = ~SA[i]] = j; } while(SA[--i] < 0);
ISAb[SA[i]] = j;
}
/* Construct the inverse suffix array of type B* suffixes using trsort. */
trsort(ISAb, SA, m, 1);
/* Set the sorted order of type B* suffixes. */
for(i = n - 1, j = m, c0 = T[n - 1]; 0 <= i;) {
for(--i, c1 = c0; (0 <= i) && ((c0 = T[i]) >= c1); --i, c1 = c0) { }
if(0 <= i) {
t = i;
for(--i, c1 = c0; (0 <= i) && ((c0 = T[i]) <= c1); --i, c1 = c0) { }
SA[ISAb[--j]] = ((t == 0) || (1 < (t - i))) ? t : ~t;
}
}
/* Calculate the index of start/end point of each bucket. */
BUCKET_B(ALPHABET_SIZE - 1, ALPHABET_SIZE - 1) = n; /* end point */
for(c0 = ALPHABET_SIZE - 2, k = m - 1; 0 <= c0; --c0) {
i = BUCKET_A(c0 + 1) - 1;
for(c1 = ALPHABET_SIZE - 1; c0 < c1; --c1) {
t = i - BUCKET_B(c0, c1);
BUCKET_B(c0, c1) = i; /* end point */
/* Move all type B* suffixes to the correct position. */
for(i = t, j = BUCKET_BSTAR(c0, c1);
j <= k;
--i, --k) { SA[i] = SA[k]; }
}
BUCKET_BSTAR(c0, c0 + 1) = i - BUCKET_B(c0, c0) + 1; /* start point */
BUCKET_B(c0, c0) = i; /* end point */
}
}
return m;
}
/* Constructs the suffix array by using the sorted order of type B* suffixes. */
static
void
construct_SA(const unsigned char *T, int *SA,
int *bucket_A, int *bucket_B,
int n, int m) {
int *i, *j, *k;
int s;
int c0, c1, c2;
if(0 < m) {
/* Construct the sorted order of type B suffixes by using
the sorted order of type B* suffixes. */
for(c1 = ALPHABET_SIZE - 2; 0 <= c1; --c1) {
/* Scan the suffix array from right to left. */
for(i = SA + BUCKET_BSTAR(c1, c1 + 1),
j = SA + BUCKET_A(c1 + 1) - 1, k = NULL, c2 = -1;
i <= j;
--j) {
if(0 < (s = *j)) {
assert(T[s] == c1);
assert(((s + 1) < n) && (T[s] <= T[s + 1]));
assert(T[s - 1] <= T[s]);
*j = ~s;
c0 = T[--s];
if((0 < s) && (T[s - 1] > c0)) { s = ~s; }
if(c0 != c2) {
if(0 <= c2) { BUCKET_B(c2, c1) = k - SA; }
k = SA + BUCKET_B(c2 = c0, c1);
}
assert(k < j); assert(k != NULL);
*k-- = s;
} else {
assert(((s == 0) && (T[s] == c1)) || (s < 0));
*j = ~s;
}
}
}
}
/* Construct the suffix array by using
the sorted order of type B suffixes. */
k = SA + BUCKET_A(c2 = T[n - 1]);
*k++ = (T[n - 2] < c2) ? ~(n - 1) : (n - 1);
/* Scan the suffix array from left to right. */
for(i = SA, j = SA + n; i < j; ++i) {
if(0 < (s = *i)) {
assert(T[s - 1] >= T[s]);
c0 = T[--s];
if((s == 0) || (T[s - 1] < c0)) { s = ~s; }
if(c0 != c2) {
BUCKET_A(c2) = k - SA;
k = SA + BUCKET_A(c2 = c0);
}
assert(i < k);
*k++ = s;
} else {
assert(s < 0);
*i = ~s;
}
}
}
/* Constructs the burrows-wheeler transformed string directly
by using the sorted order of type B* suffixes. */
static
int
construct_BWT(const unsigned char *T, int *SA,
int *bucket_A, int *bucket_B,
int n, int m) {
int *i, *j, *k, *orig;
int s;
int c0, c1, c2;
if(0 < m) {
/* Construct the sorted order of type B suffixes by using
the sorted order of type B* suffixes. */
for(c1 = ALPHABET_SIZE - 2; 0 <= c1; --c1) {
/* Scan the suffix array from right to left. */
for(i = SA + BUCKET_BSTAR(c1, c1 + 1),
j = SA + BUCKET_A(c1 + 1) - 1, k = NULL, c2 = -1;
i <= j;
--j) {
if(0 < (s = *j)) {
assert(T[s] == c1);
assert(((s + 1) < n) && (T[s] <= T[s + 1]));
assert(T[s - 1] <= T[s]);
c0 = T[--s];
*j = ~((int)c0);
if((0 < s) && (T[s - 1] > c0)) { s = ~s; }
if(c0 != c2) {
if(0 <= c2) { BUCKET_B(c2, c1) = k - SA; }
k = SA + BUCKET_B(c2 = c0, c1);
}
assert(k < j); assert(k != NULL);
*k-- = s;
} else if(s != 0) {
*j = ~s;
#ifndef NDEBUG
} else {
assert(T[s] == c1);
#endif
}
}
}
}
/* Construct the BWTed string by using
the sorted order of type B suffixes. */
k = SA + BUCKET_A(c2 = T[n - 1]);
*k++ = (T[n - 2] < c2) ? ~((int)T[n - 2]) : (n - 1);
/* Scan the suffix array from left to right. */
for(i = SA, j = SA + n, orig = SA; i < j; ++i) {
if(0 < (s = *i)) {
assert(T[s - 1] >= T[s]);
c0 = T[--s];
*i = c0;
if((0 < s) && (T[s - 1] < c0)) { s = ~((int)T[s - 1]); }
if(c0 != c2) {
BUCKET_A(c2) = k - SA;
k = SA + BUCKET_A(c2 = c0);
}
assert(i < k);
*k++ = s;
} else if(s != 0) {
*i = ~s;
} else {
orig = i;
}
}
return orig - SA;
}
/* Constructs the burrows-wheeler transformed string directly
by using the sorted order of type B* suffixes. */
static
int
construct_BWT_indexes(const unsigned char *T, int *SA,
int *bucket_A, int *bucket_B,
int n, int m,
unsigned char * num_indexes, int * indexes) {
int *i, *j, *k, *orig;
int s;
int c0, c1, c2;
int mod = n / 8;
{
mod |= mod >> 1; mod |= mod >> 2;
mod |= mod >> 4; mod |= mod >> 8;
mod |= mod >> 16; mod >>= 1;
*num_indexes = (unsigned char)((n - 1) / (mod + 1));
}
if(0 < m) {
/* Construct the sorted order of type B suffixes by using
the sorted order of type B* suffixes. */
for(c1 = ALPHABET_SIZE - 2; 0 <= c1; --c1) {
/* Scan the suffix array from right to left. */
for(i = SA + BUCKET_BSTAR(c1, c1 + 1),
j = SA + BUCKET_A(c1 + 1) - 1, k = NULL, c2 = -1;
i <= j;
--j) {
if(0 < (s = *j)) {
assert(T[s] == c1);
assert(((s + 1) < n) && (T[s] <= T[s + 1]));
assert(T[s - 1] <= T[s]);
if ((s & mod) == 0) indexes[s / (mod + 1) - 1] = j - SA;
c0 = T[--s];
*j = ~((int)c0);
if((0 < s) && (T[s - 1] > c0)) { s = ~s; }
if(c0 != c2) {
if(0 <= c2) { BUCKET_B(c2, c1) = k - SA; }
k = SA + BUCKET_B(c2 = c0, c1);
}
assert(k < j); assert(k != NULL);
*k-- = s;
} else if(s != 0) {
*j = ~s;
#ifndef NDEBUG
} else {
assert(T[s] == c1);
#endif
}
}
}
}
/* Construct the BWTed string by using
the sorted order of type B suffixes. */
k = SA + BUCKET_A(c2 = T[n - 1]);
if (T[n - 2] < c2) {
if (((n - 1) & mod) == 0) indexes[(n - 1) / (mod + 1) - 1] = k - SA;
*k++ = ~((int)T[n - 2]);
}
else {
*k++ = n - 1;
}
/* Scan the suffix array from left to right. */
for(i = SA, j = SA + n, orig = SA; i < j; ++i) {
if(0 < (s = *i)) {
assert(T[s - 1] >= T[s]);
if ((s & mod) == 0) indexes[s / (mod + 1) - 1] = i - SA;
c0 = T[--s];
*i = c0;
if(c0 != c2) {
BUCKET_A(c2) = k - SA;
k = SA + BUCKET_A(c2 = c0);
}
assert(i < k);
if((0 < s) && (T[s - 1] < c0)) {
if ((s & mod) == 0) indexes[s / (mod + 1) - 1] = k - SA;
*k++ = ~((int)T[s - 1]);
} else
*k++ = s;
} else if(s != 0) {
*i = ~s;
} else {
orig = i;
}
}
return orig - SA;
}
/*---------------------------------------------------------------------------*/
/*- Function -*/
int
divsufsort(const unsigned char *T, int *SA, int n, int openMP) {
int *bucket_A, *bucket_B;
int m;
int err = 0;
/* Check arguments. */
if((T == NULL) || (SA == NULL) || (n < 0)) { return -1; }
else if(n == 0) { return 0; }
else if(n == 1) { SA[0] = 0; return 0; }
else if(n == 2) { m = (T[0] < T[1]); SA[m ^ 1] = 0, SA[m] = 1; return 0; }
bucket_A = (int *)malloc(BUCKET_A_SIZE * sizeof(int));
bucket_B = (int *)malloc(BUCKET_B_SIZE * sizeof(int));
/* Suffixsort. */
if((bucket_A != NULL) && (bucket_B != NULL)) {
m = sort_typeBstar(T, SA, bucket_A, bucket_B, n, openMP);
construct_SA(T, SA, bucket_A, bucket_B, n, m);
} else {
err = -2;
}
free(bucket_B);
free(bucket_A);
return err;
}
int
divbwt(const unsigned char *T, unsigned char *U, int *A, int n, unsigned char * num_indexes, int * indexes, int openMP) {
int *B;
int *bucket_A, *bucket_B;
int m, pidx, i;
/* Check arguments. */
if((T == NULL) || (U == NULL) || (n < 0)) { return -1; }
else if(n <= 1) { if(n == 1) { U[0] = T[0]; } return n; }
if((B = A) == NULL) { B = (int *)malloc((size_t)(n + 1) * sizeof(int)); }
bucket_A = (int *)malloc(BUCKET_A_SIZE * sizeof(int));
bucket_B = (int *)malloc(BUCKET_B_SIZE * sizeof(int));
/* Burrows-Wheeler Transform. */
if((B != NULL) && (bucket_A != NULL) && (bucket_B != NULL)) {
m = sort_typeBstar(T, B, bucket_A, bucket_B, n, openMP);
if (num_indexes == NULL || indexes == NULL) {
pidx = construct_BWT(T, B, bucket_A, bucket_B, n, m);
} else {
pidx = construct_BWT_indexes(T, B, bucket_A, bucket_B, n, m, num_indexes, indexes);
}
/* Copy to output string. */
U[0] = T[n - 1];
for(i = 0; i < pidx; ++i) { U[i + 1] = (unsigned char)B[i]; }
for(i += 1; i < n; ++i) { U[i] = (unsigned char)B[i]; }
pidx += 1;
} else {
pidx = -2;
}
free(bucket_B);
free(bucket_A);
if(A == NULL) { free(B); }
return pidx;
}
+15 -4
View File
@@ -71,7 +71,7 @@ RUST_TARGET_32 ?= i686-unknown-linux-gnu
RUST_STATICLIB_32 := $(RUST_TARGET_DIR)/$(RUST_TARGET_32)/release/libzstd_rs.a RUST_STATICLIB_32 := $(RUST_TARGET_DIR)/$(RUST_TARGET_32)/release/libzstd_rs.a
RUST_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \ RUST_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \
--target-dir $(RUST_TARGET_DIR) --no-default-features --target-dir $(RUST_TARGET_DIR) --no-default-features
RUST_CARGO_FLAGS += --features compression,decompression RUST_CARGO_FLAGS += --features compression,decompression,dict-builder
ifneq ($(RUST_HUF_FEATURE),) ifneq ($(RUST_HUF_FEATURE),)
RUST_CARGO_FLAGS += --features $(RUST_HUF_FEATURE) RUST_CARGO_FLAGS += --features $(RUST_HUF_FEATURE)
endif endif
@@ -96,7 +96,7 @@ $(RUST_CLI_STATICLIB): $(RUST_CLI_SOURCES)
$(RUST_CLI_STATICLIB_32): $(RUST_CLI_SOURCES) $(RUST_CLI_STATICLIB_32): $(RUST_CLI_SOURCES)
$(CARGO) build $(RUST_CLI_CARGO_FLAGS) --target $(RUST_TARGET_32) $(CARGO) build $(RUST_CLI_CARGO_FLAGS) --target $(RUST_TARGET_32)
RUST_DECOMPRESS_BUILD_CONFIG := lib-c0-d1-$(RUST_BUILD_CONFIG) RUST_DECOMPRESS_BUILD_CONFIG := lib-c0-d1-b0-$(RUST_BUILD_CONFIG)
RUST_DECOMPRESS_TARGET_DIR := $(RUST_DIR)/target/$(RUST_DECOMPRESS_BUILD_CONFIG) RUST_DECOMPRESS_TARGET_DIR := $(RUST_DIR)/target/$(RUST_DECOMPRESS_BUILD_CONFIG)
RUST_DECOMPRESS_STATICLIB := $(RUST_DECOMPRESS_TARGET_DIR)/release/libzstd_rs.a RUST_DECOMPRESS_STATICLIB := $(RUST_DECOMPRESS_TARGET_DIR)/release/libzstd_rs.a
RUST_DECOMPRESS_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \ RUST_DECOMPRESS_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \
@@ -119,7 +119,7 @@ RUST_DECOMPRESS_CLI_CARGO_FLAGS := --manifest-path $(RUST_CLI_MANIFEST) --releas
$(RUST_DECOMPRESS_CLI_STATICLIB): $(RUST_CLI_SOURCES) $(RUST_DECOMPRESS_CLI_STATICLIB): $(RUST_CLI_SOURCES)
$(CARGO) build $(RUST_DECOMPRESS_CLI_CARGO_FLAGS) $(CARGO) build $(RUST_DECOMPRESS_CLI_CARGO_FLAGS)
RUST_COMPRESS_BUILD_CONFIG := lib-c1-d0-$(RUST_BUILD_CONFIG) RUST_COMPRESS_BUILD_CONFIG := lib-c1-d0-b0-$(RUST_BUILD_CONFIG)
RUST_COMPRESS_TARGET_DIR := $(RUST_DIR)/target/$(RUST_COMPRESS_BUILD_CONFIG) RUST_COMPRESS_TARGET_DIR := $(RUST_DIR)/target/$(RUST_COMPRESS_BUILD_CONFIG)
RUST_COMPRESS_STATICLIB := $(RUST_COMPRESS_TARGET_DIR)/release/libzstd_rs.a RUST_COMPRESS_STATICLIB := $(RUST_COMPRESS_TARGET_DIR)/release/libzstd_rs.a
RUST_COMPRESS_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \ RUST_COMPRESS_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \
@@ -129,6 +129,17 @@ RUST_COMPRESS_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \
$(RUST_COMPRESS_STATICLIB): $(RUST_SOURCES) $(RUST_COMPRESS_STATICLIB): $(RUST_SOURCES)
$(CARGO) build $(RUST_COMPRESS_CARGO_FLAGS) $(CARGO) build $(RUST_COMPRESS_CARGO_FLAGS)
RUST_DICTBUILDER_BUILD_CONFIG := lib-c1-d0-b1-$(RUST_BUILD_CONFIG)
RUST_DICTBUILDER_TARGET_DIR := $(RUST_DIR)/target/$(RUST_DICTBUILDER_BUILD_CONFIG)
RUST_DICTBUILDER_STATICLIB := $(RUST_DICTBUILDER_TARGET_DIR)/release/libzstd_rs.a
RUST_DICTBUILDER_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \
--target-dir $(RUST_DICTBUILDER_TARGET_DIR) \
--no-default-features \
--features compression,dict-builder
$(RUST_DICTBUILDER_STATICLIB): $(RUST_SOURCES)
$(CARGO) build $(RUST_DICTBUILDER_CARGO_FLAGS)
RUST_COMPRESS_CLI_BUILD_CONFIG := cli-c1-d0-$(RUST_BUILD_CONFIG) RUST_COMPRESS_CLI_BUILD_CONFIG := cli-c1-d0-$(RUST_BUILD_CONFIG)
RUST_COMPRESS_CLI_TARGET_DIR := $(RUST_DIR)/target/$(RUST_COMPRESS_CLI_BUILD_CONFIG) RUST_COMPRESS_CLI_TARGET_DIR := $(RUST_DIR)/target/$(RUST_COMPRESS_CLI_BUILD_CONFIG)
RUST_COMPRESS_CLI_STATICLIB := $(RUST_COMPRESS_CLI_TARGET_DIR)/release/libzstd_cli_rs.a RUST_COMPRESS_CLI_STATICLIB := $(RUST_COMPRESS_CLI_TARGET_DIR)/release/libzstd_cli_rs.a
@@ -412,7 +423,7 @@ zstd-compress: $(ZSTDLIB_COMMON_SRC) $(ZSTDLIB_COMPRESS_SRC) zstdcli.c util.c ti
## zstd-dictBuilder: executable supporting dictionary creation and compression (only) ## zstd-dictBuilder: executable supporting dictionary creation and compression (only)
CLEAN += zstd-dictBuilder CLEAN += zstd-dictBuilder
zstd-dictBuilder: $(ZSTDLIB_COMMON_SRC) $(ZSTDLIB_COMPRESS_SRC) $(ZDICT_SRC) zstdcli.c util.c timefn.c fileio.c fileio_asyncio.c dibio.c $(RUST_COMPRESS_STATICLIB) $(RUST_COMPRESS_CLI_STATICLIB) zstd-dictBuilder: $(ZSTDLIB_COMMON_SRC) $(ZSTDLIB_COMPRESS_SRC) $(ZDICT_SRC) zstdcli.c util.c timefn.c fileio.c fileio_asyncio.c dibio.c $(RUST_DICTBUILDER_STATICLIB) $(RUST_COMPRESS_CLI_STATICLIB)
$(CC) $(FLAGS) -DZSTD_NOBENCH -DZSTD_NODECOMPRESS -DZSTD_NOTRACE $^ -o $@$(EXT) $(CC) $(FLAGS) -DZSTD_NOBENCH -DZSTD_NODECOMPRESS -DZSTD_NOTRACE $^ -o $@$(EXT)
RUST_DIRECT_LINK_TARGETS := zstd32 zstd-nolegacy zstd-small zstd-frugal \ RUST_DIRECT_LINK_TARGETS := zstd32 zstd-nolegacy zstd-small zstd-frugal \
+2 -1
View File
@@ -7,9 +7,10 @@ edition = "2021"
crate-type = ["staticlib"] crate-type = ["staticlib"]
[features] [features]
default = ["compression", "decompression"] default = ["compression", "decompression", "dict-builder"]
compression = [] compression = []
decompression = [] decompression = []
dict-builder = []
huf-force-decompress-x1 = [] huf-force-decompress-x1 = []
huf-force-decompress-x2 = [] huf-force-decompress-x2 = []
+12 -6
View File
@@ -42,6 +42,10 @@ zstd ABI:
the dynamic-programming optimal parser itself remains in C for now. the dynamic-programming optimal parser itself remains in C for now.
- `zstd_ldm` implements long-distance-match parameter selection, table - `zstd_ldm` implements long-distance-match parameter selection, table
maintenance, sequence generation, and sequence consumption. maintenance, sequence generation, and sequence consumption.
- Dictionary building
- `divsufsort` constructs the suffix array that drives the legacy `ZDICT`
trainer (`ZDICT_trainFromBuffer_legacy`). The sample analysis and
dictionary assembly in `zdict.c`, `cover.c`, and `fastcover.c` remain C.
- Runtime support - Runtime support
- `threading` provides platform pthread wrappers required by zstd headers. - `threading` provides platform pthread wrappers required by zstd headers.
- `pool` implements the bounded worker pool used by multithreaded compression. - `pool` implements the bounded worker pool used by multithreaded compression.
@@ -60,10 +64,11 @@ zstd ABI:
`fileio` backend still owns file opening, safe replacement, sparse writes, `fileio` backend still owns file opening, safe replacement, sparse writes,
metadata, and streaming I/O. metadata, and streaming I/O.
The optimal block matcher, high-level frame compression, dictionary-building, The optimal block matcher, high-level frame compression, dictionary-building
legacy decoding callbacks, and the CLI file-I/O backend are still C. They must except suffix-array construction, legacy decoding callbacks, and the CLI
move before the rewrite is complete. Keeping that boundary explicit prevents a file-I/O backend are still C. They must move before the rewrite is complete.
passing hybrid build from being mistaken for the final all-Rust result. Keeping that boundary explicit prevents a passing hybrid build from being
mistaken for the final all-Rust result.
## Compatibility boundary ## Compatibility boundary
@@ -74,8 +79,9 @@ makefile source list as a small shim so header configuration and platform
preprocessor behavior stay available during the transition. preprocessor behavior stay available during the transition.
The library, test, and program makefiles select an archive directory for the The library, test, and program makefiles select an archive directory for the
active C configuration: enabled compression/decompression modules, default or active C configuration: enabled compression/decompression/dictionary-builder
forced HUF X1/X2, and the matching Rust target for 32-bit C binaries. The modules, default or forced HUF X1/X2, and the matching Rust target for 32-bit
C binaries. The
native static archive flattens Rust object members rather than nesting a Rust native static archive flattens Rust object members rather than nesting a Rust
archive, while the native shared library retains all migrated Rust exports. archive, while the native shared library retains all migrated Rust exports.
When the HUF mode changes, the test and program paths also rebuild cached C When the HUF mode changes, the test and program paths also rebuild cached C
+2756
View File
@@ -0,0 +1,2756 @@
#![allow(clippy::missing_safety_doc)]
#![allow(clippy::too_many_arguments)]
//! Suffix-array construction for the dictionary builder.
//!
//! Port of `lib/dictBuilder/divsufsort.c` (libdivsufsort-lite, Copyright (c)
//! 2003-2008 Yuta Mori, MIT license) in the exact configuration zstd compiles
//! it with: `ALPHABET_SIZE = 256`, `SS_INSERTIONSORT_THRESHOLD = 8`,
//! `SS_BLOCKSIZE = 1024`, and no OpenMP. Only `divsufsort()` is exported;
//! `divbwt()` has no callers anywhere in zstd and was not ported.
//!
//! The C implementation walks raw `int*` cursors through the caller's SA
//! buffer, including transient one-before-the-range positions. Every such
//! cursor is translated to an `isize` index into one `&mut [i32]` slice
//! covering the whole buffer, so all arithmetic — including the
//! bitwise-complement rank marking and the C `int` value semantics — matches
//! the original exactly while staying bounds-checked.
use std::os::raw::c_int;
use std::slice;
const BUCKET_A_SIZE: usize = 256; /* ALPHABET_SIZE */
const BUCKET_B_SIZE: usize = 256 * 256; /* ALPHABET_SIZE * ALPHABET_SIZE */
const ALPHABET_SIZE: i32 = 256;
const SS_INSERTIONSORT_THRESHOLD: isize = 8;
const SS_BLOCKSIZE: isize = 1024;
/* minstacksize = log(SS_BLOCKSIZE) / log(3) * 2 */
const SS_MISORT_STACKSIZE: usize = 16;
const SS_SMERGE_STACKSIZE: usize = 32;
const TR_INSERTIONSORT_THRESHOLD: isize = 8;
const TR_STACKSIZE: usize = 64;
#[rustfmt::skip]
static LG_TABLE: [i32; 256] = [
-1,0,1,1,2,2,2,2,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,
5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
];
#[rustfmt::skip]
static SQQ_TABLE: [i32; 256] = [
0, 16, 22, 27, 32, 35, 39, 42, 45, 48, 50, 53, 55, 57, 59, 61,
64, 65, 67, 69, 71, 73, 75, 76, 78, 80, 81, 83, 84, 86, 87, 89,
90, 91, 93, 94, 96, 97, 98, 99, 101, 102, 103, 104, 106, 107, 108, 109,
110, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126,
128, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142,
143, 144, 144, 145, 146, 147, 148, 149, 150, 150, 151, 152, 153, 154, 155, 155,
156, 157, 158, 159, 160, 160, 161, 162, 163, 163, 164, 165, 166, 167, 167, 168,
169, 170, 170, 171, 172, 173, 173, 174, 175, 176, 176, 177, 178, 178, 179, 180,
181, 181, 182, 183, 183, 184, 185, 185, 186, 187, 187, 188, 189, 189, 190, 191,
192, 192, 193, 193, 194, 195, 195, 196, 197, 197, 198, 199, 199, 200, 201, 201,
202, 203, 203, 204, 204, 205, 206, 206, 207, 208, 208, 209, 209, 210, 211, 211,
212, 212, 213, 214, 214, 215, 215, 216, 217, 217, 218, 218, 219, 219, 220, 221,
221, 222, 222, 223, 224, 224, 225, 225, 226, 226, 227, 227, 228, 229, 229, 230,
230, 231, 231, 232, 232, 233, 234, 234, 235, 235, 236, 236, 237, 237, 238, 238,
239, 240, 240, 241, 241, 242, 242, 243, 243, 244, 244, 245, 245, 246, 246, 247,
247, 248, 248, 249, 249, 250, 250, 251, 251, 252, 252, 253, 253, 254, 254, 255,
];
/* `ss_ilg` in its `256 <= SS_BLOCKSIZE` configuration. */
#[inline]
fn ss_ilg(n: isize) -> i32 {
let n = n as i32;
if n & 0xff00 != 0 {
8 + LG_TABLE[((n >> 8) & 0xff) as usize]
} else {
LG_TABLE[(n & 0xff) as usize]
}
}
#[inline]
fn ss_isqrt(x: isize) -> isize {
if x >= SS_BLOCKSIZE * SS_BLOCKSIZE {
return SS_BLOCKSIZE;
}
let x = x as i32;
let e = if (x as u32) & 0xffff_0000 != 0 {
if (x as u32) & 0xff00_0000 != 0 {
24 + LG_TABLE[((x >> 24) & 0xff) as usize]
} else {
16 + LG_TABLE[((x >> 16) & 0xff) as usize]
}
} else if x & 0xff00 != 0 {
8 + LG_TABLE[((x >> 8) & 0xff) as usize]
} else {
LG_TABLE[(x & 0xff) as usize]
};
let mut y;
if e >= 16 {
y = SQQ_TABLE[(x >> ((e - 6) - (e & 1))) as usize] << ((e >> 1) - 7);
if e >= 24 {
y = (y + 1 + x / y) >> 1;
}
y = (y + 1 + x / y) >> 1;
} else if e >= 8 {
y = (SQQ_TABLE[(x >> ((e - 6) - (e & 1))) as usize] >> (7 - (e >> 1))) + 1;
} else {
return (SQQ_TABLE[x as usize] >> 4) as isize;
}
(if x < y * y { y - 1 } else { y }) as isize
}
/* --------------------------------------------------------------------- */
/// Compares two suffixes. `(p10, p11)` and `(p20, p21)` are the `p[0]`/`p[1]`
/// pairs the C routine reads through its `const int*` arguments; passing the
/// values directly also serves `sssort()`'s local two-element `PAi` array.
#[inline]
fn ss_compare(t: &[u8], p10: i32, p11: i32, p20: i32, p21: i32, depth: i32) -> i32 {
let mut u1 = (depth + p10) as isize;
let mut u2 = (depth + p20) as isize;
let u1n = (p11 + 2) as isize;
let u2n = (p21 + 2) as isize;
while u1 < u1n && u2 < u2n && t[u1 as usize] == t[u2 as usize] {
u1 += 1;
u2 += 1;
}
if u1 < u1n {
if u2 < u2n {
t[u1 as usize] as i32 - t[u2 as usize] as i32
} else {
1
}
} else if u2 < u2n {
-1
} else {
0
}
}
/// `ss_compare(T, p1, p2, depth)` for pointers `p1`/`p2` into the SA buffer.
#[inline]
fn ss_compare_pa(t: &[u8], sa: &[i32], p1: isize, p2: isize, depth: i32) -> i32 {
ss_compare(
t,
sa[p1 as usize],
sa[(p1 + 1) as usize],
sa[p2 as usize],
sa[(p2 + 1) as usize],
depth,
)
}
/* --------------------------------------------------------------------- */
/* Insertionsort for small size groups */
fn ss_insertionsort(t: &[u8], sa: &mut [i32], pa: isize, first: isize, last: isize, depth: i32) {
let mut i = last - 2;
while first <= i {
let t0 = sa[i as usize];
let mut j = i + 1;
let mut r;
loop {
r = ss_compare_pa(t, sa, pa + t0 as isize, pa + sa[j as usize] as isize, depth);
if r <= 0 {
break;
}
loop {
sa[(j - 1) as usize] = sa[j as usize];
j += 1;
if !(j < last && sa[j as usize] < 0) {
break;
}
}
if last <= j {
break;
}
}
if r == 0 {
sa[j as usize] = !sa[j as usize];
}
sa[(j - 1) as usize] = t0;
i -= 1;
}
}
/* --------------------------------------------------------------------- */
/// `Td[PA[SA[p]]]` — the depth-`td` sorting key of the suffix stored at `p`.
#[inline(always)]
fn ss_key(t: &[u8], sa: &[i32], td: isize, pa: isize, p: isize) -> i32 {
t[(td + sa[(pa + sa[p as usize] as isize) as usize] as isize) as usize] as i32
}
/// `Td[v]` for an already-loaded SA element `v` (`Td[PA[v]]` in C).
#[inline(always)]
fn ss_key_of(t: &[u8], sa: &[i32], td: isize, pa: isize, v: i32) -> i32 {
t[(td + sa[(pa + v as isize) as usize] as isize) as usize] as i32
}
/// `Td[PA[SA[p]] - 1]` — the character preceding the depth-`td` key.
#[inline(always)]
fn ss_key_pred(t: &[u8], sa: &[i32], td: isize, pa: isize, p: isize) -> i32 {
t[(td + sa[(pa + sa[p as usize] as isize) as usize] as isize - 1) as usize] as i32
}
fn ss_fixdown(t: &[u8], td: isize, sa: &mut [i32], pa: isize, base: isize, i: isize, size: isize) {
let mut i = i;
let v = sa[(base + i) as usize];
let c = ss_key_of(t, sa, td, pa, v);
loop {
let mut j = 2 * i + 1;
if j >= size {
break;
}
let mut k = j;
j += 1;
let mut d = ss_key(t, sa, td, pa, base + k);
let e = ss_key(t, sa, td, pa, base + j);
if d < e {
k = j;
d = e;
}
if d <= c {
break;
}
sa[(base + i) as usize] = sa[(base + k) as usize];
i = k;
}
sa[(base + i) as usize] = v;
}
/* Simple top-down heapsort. */
fn ss_heapsort(t: &[u8], td: isize, sa: &mut [i32], pa: isize, base: isize, size: isize) {
let mut m = size;
if size % 2 == 0 {
m -= 1;
if ss_key(t, sa, td, pa, base + m / 2) < ss_key(t, sa, td, pa, base + m) {
sa.swap((base + m) as usize, (base + m / 2) as usize);
}
}
let mut i = m / 2 - 1;
while 0 <= i {
ss_fixdown(t, td, sa, pa, base, i, m);
i -= 1;
}
if size % 2 == 0 {
sa.swap(base as usize, (base + m) as usize);
ss_fixdown(t, td, sa, pa, base, 0, m);
}
let mut i = m - 1;
while 0 < i {
let t0 = sa[base as usize];
sa[base as usize] = sa[(base + i) as usize];
ss_fixdown(t, td, sa, pa, base, 0, i);
sa[(base + i) as usize] = t0;
i -= 1;
}
}
/* --------------------------------------------------------------------- */
/* Returns the median of three elements. */
#[inline]
fn ss_median3(
t: &[u8],
sa: &[i32],
td: isize,
pa: isize,
v1: isize,
v2: isize,
v3: isize,
) -> isize {
let mut v1 = v1;
let mut v2 = v2;
if ss_key(t, sa, td, pa, v1) > ss_key(t, sa, td, pa, v2) {
std::mem::swap(&mut v1, &mut v2);
}
if ss_key(t, sa, td, pa, v2) > ss_key(t, sa, td, pa, v3) {
if ss_key(t, sa, td, pa, v1) > ss_key(t, sa, td, pa, v3) {
return v1;
}
return v3;
}
v2
}
/* Returns the median of five elements. */
#[inline]
fn ss_median5(
t: &[u8],
sa: &[i32],
td: isize,
pa: isize,
v1: isize,
v2: isize,
v3: isize,
v4: isize,
v5: isize,
) -> isize {
let mut v1 = v1;
let mut v2 = v2;
let mut v3 = v3;
let mut v4 = v4;
let mut v5 = v5;
if ss_key(t, sa, td, pa, v2) > ss_key(t, sa, td, pa, v3) {
std::mem::swap(&mut v2, &mut v3);
}
if ss_key(t, sa, td, pa, v4) > ss_key(t, sa, td, pa, v5) {
std::mem::swap(&mut v4, &mut v5);
}
if ss_key(t, sa, td, pa, v2) > ss_key(t, sa, td, pa, v4) {
std::mem::swap(&mut v2, &mut v4);
std::mem::swap(&mut v3, &mut v5);
}
if ss_key(t, sa, td, pa, v1) > ss_key(t, sa, td, pa, v3) {
std::mem::swap(&mut v1, &mut v3);
}
if ss_key(t, sa, td, pa, v1) > ss_key(t, sa, td, pa, v4) {
std::mem::swap(&mut v1, &mut v4);
std::mem::swap(&mut v3, &mut v5);
}
if ss_key(t, sa, td, pa, v3) > ss_key(t, sa, td, pa, v4) {
return v4;
}
v3
}
/* Returns the pivot element. */
#[inline]
fn ss_pivot(t: &[u8], sa: &[i32], td: isize, pa: isize, first: isize, last: isize) -> isize {
let mut t0 = last - first;
let middle = first + t0 / 2;
if t0 <= 512 {
if t0 <= 32 {
return ss_median3(t, sa, td, pa, first, middle, last - 1);
}
t0 >>= 2;
return ss_median5(
t,
sa,
td,
pa,
first,
first + t0,
middle,
last - 1 - t0,
last - 1,
);
}
t0 >>= 3;
let first = ss_median3(t, sa, td, pa, first, first + t0, first + (t0 << 1));
let middle = ss_median3(t, sa, td, pa, middle - t0, middle, middle + t0);
let last = ss_median3(t, sa, td, pa, last - 1 - (t0 << 1), last - 1 - t0, last - 1);
ss_median3(t, sa, td, pa, first, middle, last)
}
/* --------------------------------------------------------------------- */
/* Binary partition for substrings. */
/* The `>= x + 1` comparison deliberately mirrors the C expression shape. */
#[allow(clippy::int_plus_one)]
fn ss_partition(sa: &mut [i32], pa: isize, first: isize, last: isize, depth: i32) -> isize {
let mut a = first - 1;
let mut b = last;
loop {
loop {
a += 1;
if !(a < b) {
break;
}
if !(sa[(pa + sa[a as usize] as isize) as usize] + depth
>= sa[(pa + sa[a as usize] as isize + 1) as usize] + 1)
{
break;
}
sa[a as usize] = !sa[a as usize];
}
loop {
b -= 1;
if !(a < b) {
break;
}
if !(sa[(pa + sa[b as usize] as isize) as usize] + depth
< sa[(pa + sa[b as usize] as isize + 1) as usize] + 1)
{
break;
}
}
if b <= a {
break;
}
let t0 = !sa[b as usize];
sa[b as usize] = sa[a as usize];
sa[a as usize] = t0;
}
if first < a {
sa[first as usize] = !sa[first as usize];
}
a
}
/* Multikey introsort for medium size groups. */
fn ss_mintrosort(t: &[u8], sa: &mut [i32], pa: isize, first: isize, last: isize, depth: i32) {
let mut stack = [(0isize, 0isize, 0i32, 0i32); SS_MISORT_STACKSIZE];
let mut ssize = 0usize;
let mut first = first;
let mut last = last;
let mut depth = depth;
let mut limit = ss_ilg(last - first);
let mut x: i32 = 0;
loop {
if last - first <= SS_INSERTIONSORT_THRESHOLD {
if 1 < last - first {
ss_insertionsort(t, sa, pa, first, last, depth);
}
/* STACK_POP */
if ssize == 0 {
return;
}
ssize -= 1;
(first, last, depth, limit) = stack[ssize];
continue;
}
let td = depth as isize;
if limit == 0 {
ss_heapsort(t, td, sa, pa, first, last - first);
}
limit -= 1;
if limit < 0 {
let mut a = first + 1;
let mut v = ss_key(t, sa, td, pa, first);
while a < last {
x = ss_key(t, sa, td, pa, a);
if x != v {
if 1 < a - first {
break;
}
v = x;
first = a;
}
a += 1;
}
if ss_key_pred(t, sa, td, pa, first) < v {
first = ss_partition(sa, pa, first, a, depth);
}
if a - first <= last - a {
if 1 < a - first {
stack[ssize] = (a, last, depth, -1);
ssize += 1;
last = a;
depth += 1;
limit = ss_ilg(a - first);
} else {
first = a;
limit = -1;
}
} else if 1 < last - a {
stack[ssize] = (first, a, depth + 1, ss_ilg(a - first));
ssize += 1;
first = a;
limit = -1;
} else {
last = a;
depth += 1;
limit = ss_ilg(a - first);
}
continue;
}
/* choose pivot */
let mut a = ss_pivot(t, sa, td, pa, first, last);
let v = ss_key(t, sa, td, pa, a);
sa.swap(first as usize, a as usize);
/* partition */
let mut b = first;
loop {
b += 1;
if !(b < last) {
break;
}
x = ss_key(t, sa, td, pa, b);
if x != v {
break;
}
}
a = b;
if a < last && x < v {
loop {
b += 1;
if !(b < last) {
break;
}
x = ss_key(t, sa, td, pa, b);
if !(x <= v) {
break;
}
if x == v {
sa.swap(b as usize, a as usize);
a += 1;
}
}
}
let mut c = last;
loop {
c -= 1;
if !(b < c) {
break;
}
x = ss_key(t, sa, td, pa, c);
if x != v {
break;
}
}
let mut d = c;
if b < d && x > v {
loop {
c -= 1;
if !(b < c) {
break;
}
x = ss_key(t, sa, td, pa, c);
if !(x >= v) {
break;
}
if x == v {
sa.swap(c as usize, d as usize);
d -= 1;
}
}
}
while b < c {
sa.swap(b as usize, c as usize);
loop {
b += 1;
if !(b < c) {
break;
}
x = ss_key(t, sa, td, pa, b);
if !(x <= v) {
break;
}
if x == v {
sa.swap(b as usize, a as usize);
a += 1;
}
}
loop {
c -= 1;
if !(b < c) {
break;
}
x = ss_key(t, sa, td, pa, c);
if !(x >= v) {
break;
}
if x == v {
sa.swap(c as usize, d as usize);
d -= 1;
}
}
}
if a <= d {
c = b - 1;
let mut s = a - first;
let t0 = b - a;
if s > t0 {
s = t0;
}
let mut e = first;
let mut f = b - s;
while 0 < s {
sa.swap(e as usize, f as usize);
s -= 1;
e += 1;
f += 1;
}
let mut s = d - c;
let t0 = last - d - 1;
if s > t0 {
s = t0;
}
let mut e = b;
let mut f = last - s;
while 0 < s {
sa.swap(e as usize, f as usize);
s -= 1;
e += 1;
f += 1;
}
a = first + (b - a);
c = last - (d - c);
b = if v <= ss_key_pred(t, sa, td, pa, a) {
a
} else {
ss_partition(sa, pa, a, c, depth)
};
if a - first <= last - c {
if last - c <= c - b {
stack[ssize] = (b, c, depth + 1, ss_ilg(c - b));
ssize += 1;
stack[ssize] = (c, last, depth, limit);
ssize += 1;
last = a;
} else if a - first <= c - b {
stack[ssize] = (c, last, depth, limit);
ssize += 1;
stack[ssize] = (b, c, depth + 1, ss_ilg(c - b));
ssize += 1;
last = a;
} else {
stack[ssize] = (c, last, depth, limit);
ssize += 1;
stack[ssize] = (first, a, depth, limit);
ssize += 1;
first = b;
last = c;
depth += 1;
limit = ss_ilg(c - b);
}
} else if a - first <= c - b {
stack[ssize] = (b, c, depth + 1, ss_ilg(c - b));
ssize += 1;
stack[ssize] = (first, a, depth, limit);
ssize += 1;
first = c;
} else if last - c <= c - b {
stack[ssize] = (first, a, depth, limit);
ssize += 1;
stack[ssize] = (b, c, depth + 1, ss_ilg(c - b));
ssize += 1;
first = c;
} else {
stack[ssize] = (first, a, depth, limit);
ssize += 1;
stack[ssize] = (c, last, depth, limit);
ssize += 1;
first = b;
last = c;
depth += 1;
limit = ss_ilg(c - b);
}
} else {
limit += 1;
if ss_key_pred(t, sa, td, pa, first) < v {
first = ss_partition(sa, pa, first, last, depth);
limit = ss_ilg(last - first);
}
depth += 1;
}
}
}
/* --------------------------------------------------------------------- */
#[inline]
fn ss_blockswap(sa: &mut [i32], a: isize, b: isize, n: isize) {
let mut a = a;
let mut b = b;
let mut n = n;
while 0 < n {
sa.swap(a as usize, b as usize);
n -= 1;
a += 1;
b += 1;
}
}
#[inline]
fn ss_rotate(sa: &mut [i32], first: isize, middle: isize, last: isize) {
let mut first = first;
let mut last = last;
let mut l = middle - first;
let mut r = last - middle;
while 0 < l && 0 < r {
if l == r {
ss_blockswap(sa, first, middle, l);
break;
}
if l < r {
let mut a = last - 1;
let mut b = middle - 1;
let mut t0 = sa[a as usize];
loop {
sa[a as usize] = sa[b as usize];
a -= 1;
sa[b as usize] = sa[a as usize];
b -= 1;
if b < first {
sa[a as usize] = t0;
last = a;
r -= l + 1;
if r <= l {
break;
}
a -= 1;
b = middle - 1;
t0 = sa[a as usize];
}
}
} else {
let mut a = first;
let mut b = middle;
let mut t0 = sa[a as usize];
loop {
sa[a as usize] = sa[b as usize];
a += 1;
sa[b as usize] = sa[a as usize];
b += 1;
if last <= b {
sa[a as usize] = t0;
first = a + 1;
l -= r + 1;
if l <= r {
break;
}
a += 1;
b = middle;
t0 = sa[a as usize];
}
}
}
}
}
/* --------------------------------------------------------------------- */
fn ss_inplacemerge(
t: &[u8],
sa: &mut [i32],
pa: isize,
first: isize,
middle: isize,
last: isize,
depth: i32,
) {
let mut middle = middle;
let mut last = last;
loop {
let x: i32;
let p: isize;
if sa[(last - 1) as usize] < 0 {
x = 1;
p = pa + (!sa[(last - 1) as usize]) as isize;
} else {
x = 0;
p = pa + sa[(last - 1) as usize] as isize;
}
let mut a = first;
let mut len = middle - first;
let mut half = len >> 1;
let mut r: i32 = -1;
while 0 < len {
let b = a + half;
let bv = sa[b as usize];
let q = ss_compare_pa(
t,
sa,
pa + (if 0 <= bv { bv } else { !bv }) as isize,
p,
depth,
);
if q < 0 {
a = b + 1;
half -= (len & 1) ^ 1;
} else {
r = q;
}
len = half;
half >>= 1;
}
if a < middle {
if r == 0 {
sa[a as usize] = !sa[a as usize];
}
ss_rotate(sa, a, middle, last);
last -= middle - a;
middle = a;
if first == middle {
break;
}
}
last -= 1;
if x != 0 {
loop {
last -= 1;
if !(sa[last as usize] < 0) {
break;
}
}
}
if middle == last {
break;
}
}
}
/* --------------------------------------------------------------------- */
/* Merge-forward with internal buffer. */
fn ss_mergeforward(
t: &[u8],
sa: &mut [i32],
pa: isize,
first: isize,
middle: isize,
last: isize,
buf: isize,
depth: i32,
) {
let bufend = buf + (middle - first) - 1;
ss_blockswap(sa, buf, first, middle - first);
let mut a = first;
let t0 = sa[a as usize];
let mut b = buf;
let mut c = middle;
loop {
let r = ss_compare_pa(
t,
sa,
pa + sa[b as usize] as isize,
pa + sa[c as usize] as isize,
depth,
);
if r < 0 {
loop {
sa[a as usize] = sa[b as usize];
a += 1;
if bufend <= b {
sa[bufend as usize] = t0;
return;
}
sa[b as usize] = sa[a as usize];
b += 1;
if !(sa[b as usize] < 0) {
break;
}
}
} else if r > 0 {
loop {
sa[a as usize] = sa[c as usize];
a += 1;
sa[c as usize] = sa[a as usize];
c += 1;
if last <= c {
while b < bufend {
sa[a as usize] = sa[b as usize];
a += 1;
sa[b as usize] = sa[a as usize];
b += 1;
}
sa[a as usize] = sa[b as usize];
sa[b as usize] = t0;
return;
}
if !(sa[c as usize] < 0) {
break;
}
}
} else {
sa[c as usize] = !sa[c as usize];
loop {
sa[a as usize] = sa[b as usize];
a += 1;
if bufend <= b {
sa[bufend as usize] = t0;
return;
}
sa[b as usize] = sa[a as usize];
b += 1;
if !(sa[b as usize] < 0) {
break;
}
}
loop {
sa[a as usize] = sa[c as usize];
a += 1;
sa[c as usize] = sa[a as usize];
c += 1;
if last <= c {
while b < bufend {
sa[a as usize] = sa[b as usize];
a += 1;
sa[b as usize] = sa[a as usize];
b += 1;
}
sa[a as usize] = sa[b as usize];
sa[b as usize] = t0;
return;
}
if !(sa[c as usize] < 0) {
break;
}
}
}
}
}
/* Merge-backward with internal buffer. */
fn ss_mergebackward(
t: &[u8],
sa: &mut [i32],
pa: isize,
first: isize,
middle: isize,
last: isize,
buf: isize,
depth: i32,
) {
let bufend = buf + (last - middle) - 1;
ss_blockswap(sa, buf, middle, last - middle);
let mut x = 0i32;
let mut p1: isize;
let mut p2: isize;
if sa[bufend as usize] < 0 {
p1 = pa + (!sa[bufend as usize]) as isize;
x |= 1;
} else {
p1 = pa + sa[bufend as usize] as isize;
}
if sa[(middle - 1) as usize] < 0 {
p2 = pa + (!sa[(middle - 1) as usize]) as isize;
x |= 2;
} else {
p2 = pa + sa[(middle - 1) as usize] as isize;
}
let mut a = last - 1;
let t0 = sa[a as usize];
let mut b = bufend;
let mut c = middle - 1;
loop {
let r = ss_compare_pa(t, sa, p1, p2, depth);
if 0 < r {
if x & 1 != 0 {
loop {
sa[a as usize] = sa[b as usize];
a -= 1;
sa[b as usize] = sa[a as usize];
b -= 1;
if !(sa[b as usize] < 0) {
break;
}
}
x ^= 1;
}
sa[a as usize] = sa[b as usize];
a -= 1;
if b <= buf {
sa[buf as usize] = t0;
break;
}
sa[b as usize] = sa[a as usize];
b -= 1;
if sa[b as usize] < 0 {
p1 = pa + (!sa[b as usize]) as isize;
x |= 1;
} else {
p1 = pa + sa[b as usize] as isize;
}
} else if r < 0 {
if x & 2 != 0 {
loop {
sa[a as usize] = sa[c as usize];
a -= 1;
sa[c as usize] = sa[a as usize];
c -= 1;
if !(sa[c as usize] < 0) {
break;
}
}
x ^= 2;
}
sa[a as usize] = sa[c as usize];
a -= 1;
sa[c as usize] = sa[a as usize];
c -= 1;
if c < first {
while buf < b {
sa[a as usize] = sa[b as usize];
a -= 1;
sa[b as usize] = sa[a as usize];
b -= 1;
}
sa[a as usize] = sa[b as usize];
sa[b as usize] = t0;
break;
}
if sa[c as usize] < 0 {
p2 = pa + (!sa[c as usize]) as isize;
x |= 2;
} else {
p2 = pa + sa[c as usize] as isize;
}
} else {
if x & 1 != 0 {
loop {
sa[a as usize] = sa[b as usize];
a -= 1;
sa[b as usize] = sa[a as usize];
b -= 1;
if !(sa[b as usize] < 0) {
break;
}
}
x ^= 1;
}
sa[a as usize] = !sa[b as usize];
a -= 1;
if b <= buf {
sa[buf as usize] = t0;
break;
}
sa[b as usize] = sa[a as usize];
b -= 1;
if x & 2 != 0 {
loop {
sa[a as usize] = sa[c as usize];
a -= 1;
sa[c as usize] = sa[a as usize];
c -= 1;
if !(sa[c as usize] < 0) {
break;
}
}
x ^= 2;
}
sa[a as usize] = sa[c as usize];
a -= 1;
sa[c as usize] = sa[a as usize];
c -= 1;
if c < first {
while buf < b {
sa[a as usize] = sa[b as usize];
a -= 1;
sa[b as usize] = sa[a as usize];
b -= 1;
}
sa[a as usize] = sa[b as usize];
sa[b as usize] = t0;
break;
}
if sa[b as usize] < 0 {
p1 = pa + (!sa[b as usize]) as isize;
x |= 1;
} else {
p1 = pa + sa[b as usize] as isize;
}
if sa[c as usize] < 0 {
p2 = pa + (!sa[c as usize]) as isize;
x |= 2;
} else {
p2 = pa + sa[c as usize] as isize;
}
}
}
}
/// `GETIDX` — undoes the "already merged" complement marking.
#[inline(always)]
fn getidx(a: i32) -> i32 {
if 0 <= a {
a
} else {
!a
}
}
/// `MERGE_CHECK` — restores or sets the complement marks after a merge.
#[inline]
fn ss_merge_check(t: &[u8], sa: &mut [i32], pa: isize, a: isize, b: isize, c: i32, depth: i32) {
if (c & 1) != 0
|| ((c & 2) != 0
&& ss_compare_pa(
t,
sa,
pa + getidx(sa[(a - 1) as usize]) as isize,
pa + sa[a as usize] as isize,
depth,
) == 0)
{
sa[a as usize] = !sa[a as usize];
}
if (c & 4) != 0
&& ss_compare_pa(
t,
sa,
pa + getidx(sa[(b - 1) as usize]) as isize,
pa + sa[b as usize] as isize,
depth,
) == 0
{
sa[b as usize] = !sa[b as usize];
}
}
/* D&C based merge. */
fn ss_swapmerge(
t: &[u8],
sa: &mut [i32],
pa: isize,
first: isize,
middle: isize,
last: isize,
buf: isize,
bufsize: isize,
depth: i32,
) {
let mut stack = [(0isize, 0isize, 0isize, 0i32); SS_SMERGE_STACKSIZE];
let mut ssize = 0usize;
let mut first = first;
let mut middle = middle;
let mut last = last;
let mut check = 0i32;
loop {
if last - middle <= bufsize {
if first < middle && middle < last {
ss_mergebackward(t, sa, pa, first, middle, last, buf, depth);
}
ss_merge_check(t, sa, pa, first, last, check, depth);
if ssize == 0 {
return;
}
ssize -= 1;
(first, middle, last, check) = stack[ssize];
continue;
}
if middle - first <= bufsize {
if first < middle {
ss_mergeforward(t, sa, pa, first, middle, last, buf, depth);
}
ss_merge_check(t, sa, pa, first, last, check, depth);
if ssize == 0 {
return;
}
ssize -= 1;
(first, middle, last, check) = stack[ssize];
continue;
}
let mut m: isize = 0;
let mut len = std::cmp::min(middle - first, last - middle);
let mut half = len >> 1;
while 0 < len {
if ss_compare_pa(
t,
sa,
pa + getidx(sa[(middle + m + half) as usize]) as isize,
pa + getidx(sa[(middle - m - half - 1) as usize]) as isize,
depth,
) < 0
{
m += half + 1;
half -= (len & 1) ^ 1;
}
len = half;
half >>= 1;
}
if 0 < m {
let lm = middle - m;
let rm = middle + m;
ss_blockswap(sa, lm, middle, m);
let mut l = middle;
let mut r = middle;
let mut next = 0i32;
if rm < last {
if sa[rm as usize] < 0 {
sa[rm as usize] = !sa[rm as usize];
if first < lm {
loop {
l -= 1;
if !(sa[l as usize] < 0) {
break;
}
}
next |= 4;
}
next |= 1;
} else if first < lm {
while sa[r as usize] < 0 {
r += 1;
}
next |= 2;
}
}
if l - first <= last - r {
stack[ssize] = (r, rm, last, (next & 3) | (check & 4));
ssize += 1;
middle = lm;
last = l;
check = (check & 3) | (next & 4);
} else {
if (next & 2) != 0 && r == middle {
next ^= 6;
}
stack[ssize] = (first, lm, l, (check & 3) | (next & 4));
ssize += 1;
first = r;
middle = rm;
check = (next & 3) | (check & 4);
}
} else {
if ss_compare_pa(
t,
sa,
pa + getidx(sa[(middle - 1) as usize]) as isize,
pa + sa[middle as usize] as isize,
depth,
) == 0
{
sa[middle as usize] = !sa[middle as usize];
}
ss_merge_check(t, sa, pa, first, last, check, depth);
if ssize == 0 {
return;
}
ssize -= 1;
(first, middle, last, check) = stack[ssize];
}
}
}
/* --------------------------------------------------------------------- */
/* Substring sort */
fn sssort(
t: &[u8],
sa: &mut [i32],
pa: isize,
first: isize,
last: isize,
buf: isize,
bufsize: isize,
depth: i32,
n: isize,
lastsuffix: bool,
) {
let mut first = first;
let mut buf = buf;
let mut bufsize = bufsize;
if lastsuffix {
first += 1;
}
let mut limit: isize = 0;
let mut middle = last;
if bufsize < SS_BLOCKSIZE && bufsize < last - first {
limit = ss_isqrt(last - first);
if bufsize < limit {
if SS_BLOCKSIZE < limit {
limit = SS_BLOCKSIZE;
}
middle = last - limit;
buf = middle;
bufsize = limit;
} else {
limit = 0;
}
}
let mut a = first;
let mut i: isize = 0;
while SS_BLOCKSIZE < middle - a {
ss_mintrosort(t, sa, pa, a, a + SS_BLOCKSIZE, depth);
let mut curbufsize = last - (a + SS_BLOCKSIZE);
let mut curbuf = a + SS_BLOCKSIZE;
if curbufsize <= bufsize {
curbufsize = bufsize;
curbuf = buf;
}
let mut b = a;
let mut k = SS_BLOCKSIZE;
let mut j = i;
while j & 1 != 0 {
ss_swapmerge(t, sa, pa, b - k, b, b + k, curbuf, curbufsize, depth);
b -= k;
k <<= 1;
j >>= 1;
}
a += SS_BLOCKSIZE;
i += 1;
}
ss_mintrosort(t, sa, pa, a, middle, depth);
let mut k = SS_BLOCKSIZE;
while i != 0 {
if i & 1 != 0 {
ss_swapmerge(t, sa, pa, a - k, a, middle, buf, bufsize, depth);
a -= k;
}
k <<= 1;
i >>= 1;
}
if limit != 0 {
ss_mintrosort(t, sa, pa, middle, last, depth);
ss_inplacemerge(t, sa, pa, first, middle, last, depth);
}
if lastsuffix {
/* Insert last type B* suffix. */
let pai0 = sa[(pa + sa[(first - 1) as usize] as isize) as usize];
let pai1 = (n - 2) as i32;
let i0 = sa[(first - 1) as usize];
let mut a = first;
while a < last {
let av = sa[a as usize];
if !(av < 0
|| 0 < ss_compare(
t,
pai0,
pai1,
sa[(pa + av as isize) as usize],
sa[(pa + av as isize + 1) as usize],
depth,
))
{
break;
}
sa[(a - 1) as usize] = av;
a += 1;
}
sa[(a - 1) as usize] = i0;
}
}
/* --------------------------------------------------------------------- */
#[inline]
fn tr_ilg(n: isize) -> i32 {
let n = n as i32;
if (n as u32) & 0xffff_0000 != 0 {
if (n as u32) & 0xff00_0000 != 0 {
24 + LG_TABLE[((n >> 24) & 0xff) as usize]
} else {
16 + LG_TABLE[((n >> 16) & 0xff) as usize]
}
} else if n & 0xff00 != 0 {
8 + LG_TABLE[((n >> 8) & 0xff) as usize]
} else {
LG_TABLE[(n & 0xff) as usize]
}
}
/* --------------------------------------------------------------------- */
/// `ISAd[SA[p]]` — the depth-offset rank of the suffix stored at `p`.
#[inline(always)]
fn tr_key(sa: &[i32], isad: isize, p: isize) -> i32 {
sa[(isad + sa[p as usize] as isize) as usize]
}
/* Simple insertionsort for small size groups. */
fn tr_insertionsort(sa: &mut [i32], isad: isize, first: isize, last: isize) {
let mut a = first + 1;
while a < last {
let t0 = sa[a as usize];
let mut b = a - 1;
let mut r;
loop {
r = sa[(isad + t0 as isize) as usize] - tr_key(sa, isad, b);
if !(0 > r) {
break;
}
loop {
sa[(b + 1) as usize] = sa[b as usize];
b -= 1;
if !(first <= b && sa[b as usize] < 0) {
break;
}
}
if b < first {
break;
}
}
if r == 0 {
sa[b as usize] = !sa[b as usize];
}
sa[(b + 1) as usize] = t0;
a += 1;
}
}
/* --------------------------------------------------------------------- */
fn tr_fixdown(sa: &mut [i32], isad: isize, base: isize, i: isize, size: isize) {
let mut i = i;
let v = sa[(base + i) as usize];
let c = sa[(isad + v as isize) as usize];
loop {
let mut j = 2 * i + 1;
if j >= size {
break;
}
let mut k = j;
j += 1;
let mut d = tr_key(sa, isad, base + k);
let e = tr_key(sa, isad, base + j);
if d < e {
k = j;
d = e;
}
if d <= c {
break;
}
sa[(base + i) as usize] = sa[(base + k) as usize];
i = k;
}
sa[(base + i) as usize] = v;
}
/* Simple top-down heapsort. */
fn tr_heapsort(sa: &mut [i32], isad: isize, base: isize, size: isize) {
let mut m = size;
if size % 2 == 0 {
m -= 1;
if tr_key(sa, isad, base + m / 2) < tr_key(sa, isad, base + m) {
sa.swap((base + m) as usize, (base + m / 2) as usize);
}
}
let mut i = m / 2 - 1;
while 0 <= i {
tr_fixdown(sa, isad, base, i, m);
i -= 1;
}
if size % 2 == 0 {
sa.swap(base as usize, (base + m) as usize);
tr_fixdown(sa, isad, base, 0, m);
}
let mut i = m - 1;
while 0 < i {
let t0 = sa[base as usize];
sa[base as usize] = sa[(base + i) as usize];
tr_fixdown(sa, isad, base, 0, i);
sa[(base + i) as usize] = t0;
i -= 1;
}
}
/* --------------------------------------------------------------------- */
/* Returns the median of three elements. */
#[inline]
fn tr_median3(sa: &[i32], isad: isize, v1: isize, v2: isize, v3: isize) -> isize {
let mut v1 = v1;
let mut v2 = v2;
if tr_key(sa, isad, v1) > tr_key(sa, isad, v2) {
std::mem::swap(&mut v1, &mut v2);
}
if tr_key(sa, isad, v2) > tr_key(sa, isad, v3) {
if tr_key(sa, isad, v1) > tr_key(sa, isad, v3) {
return v1;
}
return v3;
}
v2
}
/* Returns the median of five elements. */
#[inline]
fn tr_median5(
sa: &[i32],
isad: isize,
v1: isize,
v2: isize,
v3: isize,
v4: isize,
v5: isize,
) -> isize {
let mut v1 = v1;
let mut v2 = v2;
let mut v3 = v3;
let mut v4 = v4;
let mut v5 = v5;
if tr_key(sa, isad, v2) > tr_key(sa, isad, v3) {
std::mem::swap(&mut v2, &mut v3);
}
if tr_key(sa, isad, v4) > tr_key(sa, isad, v5) {
std::mem::swap(&mut v4, &mut v5);
}
if tr_key(sa, isad, v2) > tr_key(sa, isad, v4) {
std::mem::swap(&mut v2, &mut v4);
std::mem::swap(&mut v3, &mut v5);
}
if tr_key(sa, isad, v1) > tr_key(sa, isad, v3) {
std::mem::swap(&mut v1, &mut v3);
}
if tr_key(sa, isad, v1) > tr_key(sa, isad, v4) {
std::mem::swap(&mut v1, &mut v4);
std::mem::swap(&mut v3, &mut v5);
}
if tr_key(sa, isad, v3) > tr_key(sa, isad, v4) {
return v4;
}
v3
}
/* Returns the pivot element. */
#[inline]
fn tr_pivot(sa: &[i32], isad: isize, first: isize, last: isize) -> isize {
let mut t0 = last - first;
let middle = first + t0 / 2;
if t0 <= 512 {
if t0 <= 32 {
return tr_median3(sa, isad, first, middle, last - 1);
}
t0 >>= 2;
return tr_median5(sa, isad, first, first + t0, middle, last - 1 - t0, last - 1);
}
t0 >>= 3;
let first = tr_median3(sa, isad, first, first + t0, first + (t0 << 1));
let middle = tr_median3(sa, isad, middle - t0, middle, middle + t0);
let last = tr_median3(sa, isad, last - 1 - (t0 << 1), last - 1 - t0, last - 1);
tr_median3(sa, isad, first, middle, last)
}
/* --------------------------------------------------------------------- */
struct TrBudget {
chance: i32,
remain: i32,
incval: i32,
count: i32,
}
impl TrBudget {
fn new(chance: i32, incval: i32) -> Self {
TrBudget {
chance,
remain: incval,
incval,
count: 0,
}
}
fn check(&mut self, size: isize) -> bool {
let size = size as i32;
if size <= self.remain {
self.remain -= size;
return true;
}
if self.chance == 0 {
self.count += size;
return false;
}
self.remain += self.incval - size;
self.chance -= 1;
true
}
}
/* --------------------------------------------------------------------- */
fn tr_partition(
sa: &mut [i32],
isad: isize,
first: isize,
middle: isize,
last: isize,
v: i32,
) -> (isize, isize) {
let mut first = first;
let mut last = last;
let mut x: i32 = 0;
let mut b = middle - 1;
loop {
b += 1;
if !(b < last) {
break;
}
x = tr_key(sa, isad, b);
if x != v {
break;
}
}
let mut a = b;
if a < last && x < v {
loop {
b += 1;
if !(b < last) {
break;
}
x = tr_key(sa, isad, b);
if !(x <= v) {
break;
}
if x == v {
sa.swap(b as usize, a as usize);
a += 1;
}
}
}
let mut c = last;
loop {
c -= 1;
if !(b < c) {
break;
}
x = tr_key(sa, isad, c);
if x != v {
break;
}
}
let mut d = c;
if b < d && x > v {
loop {
c -= 1;
if !(b < c) {
break;
}
x = tr_key(sa, isad, c);
if !(x >= v) {
break;
}
if x == v {
sa.swap(c as usize, d as usize);
d -= 1;
}
}
}
while b < c {
sa.swap(b as usize, c as usize);
loop {
b += 1;
if !(b < c) {
break;
}
x = tr_key(sa, isad, b);
if !(x <= v) {
break;
}
if x == v {
sa.swap(b as usize, a as usize);
a += 1;
}
}
loop {
c -= 1;
if !(b < c) {
break;
}
x = tr_key(sa, isad, c);
if !(x >= v) {
break;
}
if x == v {
sa.swap(c as usize, d as usize);
d -= 1;
}
}
}
if a <= d {
c = b - 1;
let mut s = a - first;
let t0 = b - a;
if s > t0 {
s = t0;
}
let mut e = first;
let mut f = b - s;
while 0 < s {
sa.swap(e as usize, f as usize);
s -= 1;
e += 1;
f += 1;
}
let mut s = d - c;
let t0 = last - d - 1;
if s > t0 {
s = t0;
}
let mut e = b;
let mut f = last - s;
while 0 < s {
sa.swap(e as usize, f as usize);
s -= 1;
e += 1;
f += 1;
}
first += b - a;
last -= d - c;
}
(first, last)
}
/* sort suffixes of middle partition by using sorted order of suffixes of
* left and right partition. */
fn tr_copy(
sa: &mut [i32],
isa: isize,
first: isize,
a: isize,
b: isize,
last: isize,
depth: isize,
) {
/* All cursor arithmetic is relative to the slice start, which is the C
* routine's `SA` pointer, so `x - SA` becomes plain `x`. */
let v = (b - 1) as i32;
let mut c = first;
let mut d = a - 1;
while c <= d {
let s = sa[c as usize] - depth as i32;
if 0 <= s && sa[(isa + s as isize) as usize] == v {
d += 1;
sa[d as usize] = s;
sa[(isa + s as isize) as usize] = d as i32;
}
c += 1;
}
let mut c = last - 1;
let e = d + 1;
let mut d = b;
while e < d {
let s = sa[c as usize] - depth as i32;
if 0 <= s && sa[(isa + s as isize) as usize] == v {
d -= 1;
sa[d as usize] = s;
sa[(isa + s as isize) as usize] = d as i32;
}
c -= 1;
}
}
fn tr_partialcopy(
sa: &mut [i32],
isa: isize,
first: isize,
a: isize,
b: isize,
last: isize,
depth: isize,
) {
let v = (b - 1) as i32;
let mut newrank: i32 = -1;
let mut lastrank: i32 = -1;
let mut c = first;
let mut d = a - 1;
while c <= d {
let s = sa[c as usize] - depth as i32;
if 0 <= s && sa[(isa + s as isize) as usize] == v {
d += 1;
sa[d as usize] = s;
let rank = sa[(isa + s as isize + depth) as usize];
if lastrank != rank {
lastrank = rank;
newrank = d as i32;
}
sa[(isa + s as isize) as usize] = newrank;
}
c += 1;
}
let mut lastrank: i32 = -1;
let mut e = d;
while first <= e {
let rank = sa[(isa + sa[e as usize] as isize) as usize];
if lastrank != rank {
lastrank = rank;
newrank = e as i32;
}
if newrank != rank {
sa[(isa + sa[e as usize] as isize) as usize] = newrank;
}
e -= 1;
}
let mut lastrank: i32 = -1;
let mut c = last - 1;
let e = d + 1;
let mut d = b;
while e < d {
let s = sa[c as usize] - depth as i32;
if 0 <= s && sa[(isa + s as isize) as usize] == v {
d -= 1;
sa[d as usize] = s;
let rank = sa[(isa + s as isize + depth) as usize];
if lastrank != rank {
lastrank = rank;
newrank = d as i32;
}
sa[(isa + s as isize) as usize] = newrank;
}
c -= 1;
}
}
fn tr_introsort(
sa: &mut [i32],
isa: isize,
isad: isize,
first: isize,
last: isize,
budget: &mut TrBudget,
) {
/* Stack frames are (ISAd, first, last, limit, trlink); the tandem-repeat
* copy frame stores its `(a, b)` pair in the pointer fields with a zero
* placeholder where C pushes a NULL ISAd. */
let mut stack = [(0isize, 0isize, 0isize, 0i32, 0i32); TR_STACKSIZE];
let mut ssize = 0usize;
let mut trlink: i32 = -1;
let mut isad = isad;
let mut first = first;
let mut last = last;
let incr = isad - isa;
let mut limit = tr_ilg(last - first);
loop {
if limit < 0 {
if limit == -1 {
/* tandem repeat partition */
let (a, b) = tr_partition(sa, isad - incr, first, first, last, (last - 1) as i32);
/* update ranks */
if a < last {
let v = (a - 1) as i32;
let mut c = first;
while c < a {
sa[(isa + sa[c as usize] as isize) as usize] = v;
c += 1;
}
}
if b < last {
let v = (b - 1) as i32;
let mut c = a;
while c < b {
sa[(isa + sa[c as usize] as isize) as usize] = v;
c += 1;
}
}
/* push */
if 1 < b - a {
stack[ssize] = (0, a, b, 0, 0);
ssize += 1;
stack[ssize] = (isad - incr, first, last, -2, trlink);
ssize += 1;
trlink = ssize as i32 - 2;
}
if a - first <= last - b {
if 1 < a - first {
stack[ssize] = (isad, b, last, tr_ilg(last - b), trlink);
ssize += 1;
last = a;
limit = tr_ilg(a - first);
} else if 1 < last - b {
first = b;
limit = tr_ilg(last - b);
} else {
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
} else if 1 < last - b {
stack[ssize] = (isad, first, a, tr_ilg(a - first), trlink);
ssize += 1;
first = b;
limit = tr_ilg(last - b);
} else if 1 < a - first {
last = a;
limit = tr_ilg(a - first);
} else {
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
} else if limit == -2 {
/* tandem repeat copy */
ssize -= 1;
let a = stack[ssize].1;
let b = stack[ssize].2;
if stack[ssize].3 == 0 {
tr_copy(sa, isa, first, a, b, last, isad - isa);
} else {
if 0 <= trlink {
stack[trlink as usize].3 = -1;
}
tr_partialcopy(sa, isa, first, a, b, last, isad - isa);
}
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
} else {
/* sorted partition */
if 0 <= sa[first as usize] {
let mut a = first;
loop {
sa[(isa + sa[a as usize] as isize) as usize] = a as i32;
a += 1;
if !(a < last && 0 <= sa[a as usize]) {
break;
}
}
first = a;
}
if first < last {
let mut a = first;
loop {
sa[a as usize] = !sa[a as usize];
a += 1;
if !(sa[a as usize] < 0) {
break;
}
}
let next =
if sa[(isa + sa[a as usize] as isize) as usize] != tr_key(sa, isad, a) {
tr_ilg(a - first + 1)
} else {
-1
};
a += 1;
if a < last {
let v = (a - 1) as i32;
let mut b = first;
while b < a {
sa[(isa + sa[b as usize] as isize) as usize] = v;
b += 1;
}
}
/* push */
if budget.check(a - first) {
if a - first <= last - a {
stack[ssize] = (isad, a, last, -3, trlink);
ssize += 1;
isad += incr;
last = a;
limit = next;
} else if 1 < last - a {
stack[ssize] = (isad + incr, first, a, next, trlink);
ssize += 1;
first = a;
limit = -3;
} else {
isad += incr;
last = a;
limit = next;
}
} else {
if 0 <= trlink {
stack[trlink as usize].3 = -1;
}
if 1 < last - a {
first = a;
limit = -3;
} else {
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
}
} else {
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
}
continue;
}
if last - first <= TR_INSERTIONSORT_THRESHOLD {
tr_insertionsort(sa, isad, first, last);
limit = -3;
continue;
}
/* C decrements `limit` here (`limit-- == 0`); the decrement is
* observable only on the not-taken path because the taken path
* overwrites `limit` with -3. */
if limit == 0 {
tr_heapsort(sa, isad, first, last - first);
let mut a = last - 1;
while first < a {
let x = tr_key(sa, isad, a);
let mut b = a - 1;
while first <= b && tr_key(sa, isad, b) == x {
sa[b as usize] = !sa[b as usize];
b -= 1;
}
a = b;
}
limit = -3;
continue;
}
limit -= 1;
/* choose pivot */
let a = tr_pivot(sa, isad, first, last);
sa.swap(first as usize, a as usize);
let v = tr_key(sa, isad, first);
/* partition */
let (a, b) = tr_partition(sa, isad, first, first + 1, last, v);
if last - first != b - a {
let next = if sa[(isa + sa[a as usize] as isize) as usize] != v {
tr_ilg(b - a)
} else {
-1
};
/* update ranks */
{
let vv = (a - 1) as i32;
let mut c = first;
while c < a {
sa[(isa + sa[c as usize] as isize) as usize] = vv;
c += 1;
}
}
if b < last {
let vv = (b - 1) as i32;
let mut c = a;
while c < b {
sa[(isa + sa[c as usize] as isize) as usize] = vv;
c += 1;
}
}
/* push */
if 1 < b - a && budget.check(b - a) {
if a - first <= last - b {
if last - b <= b - a {
if 1 < a - first {
stack[ssize] = (isad + incr, a, b, next, trlink);
ssize += 1;
stack[ssize] = (isad, b, last, limit, trlink);
ssize += 1;
last = a;
} else if 1 < last - b {
stack[ssize] = (isad + incr, a, b, next, trlink);
ssize += 1;
first = b;
} else {
isad += incr;
first = a;
last = b;
limit = next;
}
} else if a - first <= b - a {
if 1 < a - first {
stack[ssize] = (isad, b, last, limit, trlink);
ssize += 1;
stack[ssize] = (isad + incr, a, b, next, trlink);
ssize += 1;
last = a;
} else {
stack[ssize] = (isad, b, last, limit, trlink);
ssize += 1;
isad += incr;
first = a;
last = b;
limit = next;
}
} else {
stack[ssize] = (isad, b, last, limit, trlink);
ssize += 1;
stack[ssize] = (isad, first, a, limit, trlink);
ssize += 1;
isad += incr;
first = a;
last = b;
limit = next;
}
} else if a - first <= b - a {
if 1 < last - b {
stack[ssize] = (isad + incr, a, b, next, trlink);
ssize += 1;
stack[ssize] = (isad, first, a, limit, trlink);
ssize += 1;
first = b;
} else if 1 < a - first {
stack[ssize] = (isad + incr, a, b, next, trlink);
ssize += 1;
last = a;
} else {
isad += incr;
first = a;
last = b;
limit = next;
}
} else if last - b <= b - a {
if 1 < last - b {
stack[ssize] = (isad, first, a, limit, trlink);
ssize += 1;
stack[ssize] = (isad + incr, a, b, next, trlink);
ssize += 1;
first = b;
} else {
stack[ssize] = (isad, first, a, limit, trlink);
ssize += 1;
isad += incr;
first = a;
last = b;
limit = next;
}
} else {
stack[ssize] = (isad, first, a, limit, trlink);
ssize += 1;
stack[ssize] = (isad, b, last, limit, trlink);
ssize += 1;
isad += incr;
first = a;
last = b;
limit = next;
}
} else {
if 1 < b - a && 0 <= trlink {
stack[trlink as usize].3 = -1;
}
if a - first <= last - b {
if 1 < a - first {
stack[ssize] = (isad, b, last, limit, trlink);
ssize += 1;
last = a;
} else if 1 < last - b {
first = b;
} else {
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
} else if 1 < last - b {
stack[ssize] = (isad, first, a, limit, trlink);
ssize += 1;
first = b;
} else if 1 < a - first {
last = a;
} else {
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
}
} else if budget.check(last - first) {
limit = tr_ilg(last - first);
isad += incr;
} else {
if 0 <= trlink {
stack[trlink as usize].3 = -1;
}
if ssize == 0 {
return;
}
ssize -= 1;
(isad, first, last, limit, trlink) = stack[ssize];
}
}
}
/* --------------------------------------------------------------------- */
/* Tandem repeat sort */
fn trsort(sa: &mut [i32], isa: isize, n: isize, depth: isize) {
let mut budget = TrBudget::new(tr_ilg(n) * 2 / 3, n as i32);
/* trbudget_init(&budget, tr_ilg(n) * 3 / 4, n); */
let mut isad = isa + depth;
while -(n as i32) < sa[0] {
let mut first: isize = 0;
let mut skip: isize = 0;
let mut unsorted: i32 = 0;
loop {
let t0 = sa[first as usize];
if t0 < 0 {
first -= t0 as isize;
skip += t0 as isize;
} else {
if skip != 0 {
sa[(first + skip) as usize] = skip as i32;
skip = 0;
}
let last = sa[(isa + t0 as isize) as usize] as isize + 1;
if 1 < last - first {
budget.count = 0;
tr_introsort(sa, isa, isad, first, last, &mut budget);
if budget.count != 0 {
unsorted += budget.count;
} else {
skip = first - last;
}
} else if last - first == 1 {
skip = -1;
}
first = last;
}
if !(first < n) {
break;
}
}
if skip != 0 {
sa[(first + skip) as usize] = skip as i32;
}
if unsorted == 0 {
break;
}
isad += isad - isa;
}
}
/* --------------------------------------------------------------------- */
/// `BUCKET_B(c0, c1)` for the 256-symbol alphabet.
#[inline(always)]
fn bb(c0: i32, c1: i32) -> usize {
(((c1 as u32) << 8) | c0 as u32) as usize
}
/// `BUCKET_BSTAR(c0, c1)` for the 256-symbol alphabet.
#[inline(always)]
fn bstar(c0: i32, c1: i32) -> usize {
(((c0 as u32) << 8) | c1 as u32) as usize
}
/* Sorts suffixes of type B*. */
fn sort_type_bstar(
t: &[u8],
sa: &mut [i32],
bucket_a: &mut [i32],
bucket_b: &mut [i32],
n: isize,
) -> isize {
/* Initialize bucket arrays. */
for slot in bucket_a.iter_mut() {
*slot = 0;
}
for slot in bucket_b.iter_mut() {
*slot = 0;
}
/* Count the number of occurrences of the first one or two characters of
each type A, B and B* suffix. Moreover, store the beginning position of
all type B* suffixes into the array SA. */
let mut i = n - 1;
let mut m = n;
let mut c0 = t[(n - 1) as usize] as i32;
let mut c1;
while 0 <= i {
/* type A suffix. */
loop {
c1 = c0;
bucket_a[c1 as usize] += 1;
i -= 1;
if 0 <= i {
c0 = t[i as usize] as i32;
if c0 >= c1 {
continue;
}
}
break;
}
if 0 <= i {
/* type B* suffix. */
bucket_b[bstar(c0, c1)] += 1;
m -= 1;
sa[m as usize] = i as i32;
/* type B suffix. */
i -= 1;
c1 = c0;
while 0 <= i {
c0 = t[i as usize] as i32;
if !(c0 <= c1) {
break;
}
bucket_b[bb(c0, c1)] += 1;
i -= 1;
c1 = c0;
}
}
}
let m = n - m;
/*
note:
A type B* suffix is lexicographically smaller than a type B suffix that
begins with the same first two characters.
*/
/* Calculate the index of start/end point of each bucket. */
{
let mut i: i32 = 0;
let mut j: i32 = 0;
for c0 in 0..ALPHABET_SIZE {
let t0 = i + bucket_a[c0 as usize];
bucket_a[c0 as usize] = i + j; /* start point */
i = t0 + bucket_b[bb(c0, c0)];
for c1 in (c0 + 1)..ALPHABET_SIZE {
j += bucket_b[bstar(c0, c1)];
bucket_b[bstar(c0, c1)] = j; /* end point */
i += bucket_b[bb(c0, c1)];
}
}
}
if 0 < m {
/* Sort the type B* suffixes by their first two characters. */
let pab = n - m;
let isab = m;
let mut i = m - 2;
while 0 <= i {
let t0 = sa[(pab + i) as usize];
let c0 = t[t0 as usize] as i32;
let c1 = t[(t0 + 1) as usize] as i32;
bucket_b[bstar(c0, c1)] -= 1;
sa[bucket_b[bstar(c0, c1)] as usize] = i as i32;
i -= 1;
}
{
let t0 = sa[(pab + m - 1) as usize];
let c0 = t[t0 as usize] as i32;
let c1 = t[(t0 + 1) as usize] as i32;
bucket_b[bstar(c0, c1)] -= 1;
sa[bucket_b[bstar(c0, c1)] as usize] = (m - 1) as i32;
}
/* Sort the type B* substrings using sssort. */
let buf = m;
let bufsize = n - 2 * m;
let mut c0 = ALPHABET_SIZE - 2;
let mut j = m;
while 0 < j {
let mut c1 = ALPHABET_SIZE - 1;
while c0 < c1 {
let i = bucket_b[bstar(c0, c1)] as isize;
if 1 < j - i {
sssort(
t,
sa,
pab,
i,
j,
buf,
bufsize,
2,
n,
sa[i as usize] == (m - 1) as i32,
);
}
j = i;
c1 -= 1;
}
c0 -= 1;
}
/* Compute ranks of type B* substrings. */
let mut i = m - 1;
while 0 <= i {
if 0 <= sa[i as usize] {
let j = i;
loop {
sa[(isab + sa[i as usize] as isize) as usize] = i as i32;
i -= 1;
if !(0 <= i && 0 <= sa[i as usize]) {
break;
}
}
sa[(i + 1) as usize] = (i - j) as i32;
if i <= 0 {
break;
}
}
let j = i;
loop {
sa[i as usize] = !sa[i as usize];
sa[(isab + sa[i as usize] as isize) as usize] = j as i32;
i -= 1;
if !(sa[i as usize] < 0) {
break;
}
}
sa[(isab + sa[i as usize] as isize) as usize] = j as i32;
i -= 1;
}
/* Construct the inverse suffix array of type B* suffixes using
trsort. */
trsort(sa, isab, m, 1);
/* Set the sorted order of type B* suffixes. */
let mut i = n - 1;
let mut j = m;
let mut c0 = t[(n - 1) as usize] as i32;
while 0 <= i {
i -= 1;
let mut c1 = c0;
while 0 <= i {
c0 = t[i as usize] as i32;
if !(c0 >= c1) {
break;
}
i -= 1;
c1 = c0;
}
if 0 <= i {
let t0 = i;
i -= 1;
c1 = c0;
while 0 <= i {
c0 = t[i as usize] as i32;
if !(c0 <= c1) {
break;
}
i -= 1;
c1 = c0;
}
j -= 1;
sa[sa[(isab + j) as usize] as usize] = if t0 == 0 || 1 < t0 - i {
t0 as i32
} else {
!(t0 as i32)
};
}
}
/* Calculate the index of start/end point of each bucket. */
bucket_b[bb(ALPHABET_SIZE - 1, ALPHABET_SIZE - 1)] = n as i32; /* end point */
let mut k = m - 1;
let mut c0 = ALPHABET_SIZE - 2;
while 0 <= c0 {
let mut i = bucket_a[(c0 + 1) as usize] as isize - 1;
let mut c1 = ALPHABET_SIZE - 1;
while c0 < c1 {
let t0 = i - bucket_b[bb(c0, c1)] as isize;
bucket_b[bb(c0, c1)] = i as i32; /* end point */
/* Move all type B* suffixes to the correct position. */
i = t0;
let j = bucket_b[bstar(c0, c1)] as isize;
while j <= k {
sa[i as usize] = sa[k as usize];
i -= 1;
k -= 1;
}
c1 -= 1;
}
bucket_b[bstar(c0, c0 + 1)] = (i - bucket_b[bb(c0, c0)] as isize + 1) as i32; /* start point */
bucket_b[bb(c0, c0)] = i as i32; /* end point */
c0 -= 1;
}
}
m
}
/* Constructs the suffix array by using the sorted order of type B*
* suffixes. */
fn construct_sa(
t: &[u8],
sa: &mut [i32],
bucket_a: &mut [i32],
bucket_b: &mut [i32],
n: isize,
m: isize,
) {
if 0 < m {
/* Construct the sorted order of type B suffixes by using
the sorted order of type B* suffixes. */
let mut c1 = ALPHABET_SIZE - 2;
while 0 <= c1 {
/* Scan the suffix array from right to left. */
let i = bucket_b[bstar(c1, c1 + 1)] as isize;
let mut j = bucket_a[(c1 + 1) as usize] as isize - 1;
let mut k: isize = 0;
let mut c2: i32 = -1;
while i <= j {
let mut s = sa[j as usize];
if 0 < s {
debug_assert_eq!(t[s as usize] as i32, c1);
debug_assert!((s as isize + 1) < n && t[s as usize] <= t[(s + 1) as usize]);
debug_assert!(t[(s - 1) as usize] <= t[s as usize]);
sa[j as usize] = !s;
s -= 1;
let c0 = t[s as usize] as i32;
if 0 < s && (t[(s - 1) as usize] as i32) > c0 {
s = !s;
}
if c0 != c2 {
if 0 <= c2 {
bucket_b[bb(c2, c1)] = k as i32;
}
c2 = c0;
k = bucket_b[bb(c2, c1)] as isize;
}
debug_assert!(k < j);
sa[k as usize] = s;
k -= 1;
} else {
debug_assert!((s == 0 && t[s as usize] as i32 == c1) || s < 0);
sa[j as usize] = !s;
}
j -= 1;
}
c1 -= 1;
}
}
/* Construct the suffix array by using the sorted order of type B
suffixes. */
let mut c2 = t[(n - 1) as usize] as i32;
let mut k = bucket_a[c2 as usize] as isize;
sa[k as usize] = if (t[(n - 2) as usize] as i32) < c2 {
!((n - 1) as i32)
} else {
(n - 1) as i32
};
k += 1;
/* Scan the suffix array from left to right. */
let mut i: isize = 0;
let j = n;
while i < j {
let mut s = sa[i as usize];
if 0 < s {
debug_assert!(t[(s - 1) as usize] >= t[s as usize]);
s -= 1;
let c0 = t[s as usize] as i32;
if s == 0 || (t[(s - 1) as usize] as i32) < c0 {
s = !s;
}
if c0 != c2 {
bucket_a[c2 as usize] = k as i32;
c2 = c0;
k = bucket_a[c2 as usize] as isize;
}
debug_assert!(i < k);
sa[k as usize] = s;
k += 1;
} else {
debug_assert!(s < 0);
sa[i as usize] = !s;
}
i += 1;
}
}
/* --------------------------------------------------------------------- */
/// Rust implementation of the `divsufsort()` entry point used by
/// `ZDICT_trainFromBuffer_legacy()`.
///
/// Integration removes the C function body, so this direct export provides
/// the existing library symbol without a wrapper. The `open_mp` parameter is
/// accepted for signature compatibility only: zstd never defines
/// `LIBBSC_OPENMP`, so the C implementation ignored it as well.
///
/// Returns 0 on success, -1 for invalid arguments, and -2 when the bucket
/// work arrays cannot be allocated, exactly like the C routine.
#[no_mangle]
pub unsafe extern "C" fn divsufsort(
t: *const u8,
sa: *mut c_int,
n: c_int,
open_mp: c_int,
) -> c_int {
let _ = open_mp;
/* Check arguments. */
if t.is_null() || sa.is_null() || n < 0 {
return -1;
}
if n == 0 {
return 0;
}
let text = unsafe { slice::from_raw_parts(t, n as usize) };
let suffix = unsafe { slice::from_raw_parts_mut(sa, n as usize) };
if n == 1 {
suffix[0] = 0;
return 0;
}
if n == 2 {
let m = usize::from(text[0] < text[1]);
suffix[m ^ 1] = 0;
suffix[m] = 1;
return 0;
}
let mut bucket_a: Vec<i32> = Vec::new();
let mut bucket_b: Vec<i32> = Vec::new();
if bucket_a.try_reserve_exact(BUCKET_A_SIZE).is_err()
|| bucket_b.try_reserve_exact(BUCKET_B_SIZE).is_err()
{
/* Match the C implementation's -2 result when malloc fails. */
return -2;
}
bucket_a.resize(BUCKET_A_SIZE, 0);
bucket_b.resize(BUCKET_B_SIZE, 0);
/* Suffixsort. */
let m = sort_type_bstar(text, suffix, &mut bucket_a, &mut bucket_b, n as isize);
construct_sa(text, suffix, &mut bucket_a, &mut bucket_b, n as isize, m);
0
}
#[cfg(test)]
mod tests {
use super::*;
use std::ptr;
fn build_sa(text: &[u8]) -> Vec<i32> {
let mut sa = vec![0i32; text.len()];
let result = unsafe { divsufsort(text.as_ptr(), sa.as_mut_ptr(), text.len() as c_int, 0) };
assert_eq!(result, 0);
sa
}
/// Trivial O(n^2 log n) reference: sort the suffix start positions by the
/// suffixes themselves.
fn reference_sa(text: &[u8]) -> Vec<i32> {
let mut sa: Vec<i32> = (0..text.len() as i32).collect();
sa.sort_by(|&a, &b| text[a as usize..].cmp(&text[b as usize..]));
sa
}
/// Suffix-array invariants: a permutation of `0..n` whose suffixes are in
/// strictly increasing lexicographic order.
fn assert_valid_sa(text: &[u8], sa: &[i32]) {
assert_eq!(sa.len(), text.len());
let mut seen = vec![false; text.len()];
for &p in sa {
let p = usize::try_from(p).expect("suffix index must be non-negative");
assert!(p < text.len(), "suffix index {p} out of range");
assert!(!seen[p], "duplicate suffix index {p}");
seen[p] = true;
}
for pair in sa.windows(2) {
assert!(
text[pair[0] as usize..] < text[pair[1] as usize..],
"suffixes {} and {} are not in sorted order",
pair[0],
pair[1]
);
}
}
/// Fixed-seed numerical-recipes LCG, used to generate reproducible
/// pseudo-random sample buffers.
fn lcg_bytes(len: usize, seed: u32, alphabet: u32) -> Vec<u8> {
let mut state = seed;
(0..len)
.map(|_| {
state = state.wrapping_mul(1_664_525).wrapping_add(1_013_904_223);
((state >> 24) % alphabet) as u8
})
.collect()
}
#[test]
fn rejects_invalid_arguments() {
let text = [0u8; 1];
let mut sa = [0i32; 1];
assert_eq!(
unsafe { divsufsort(ptr::null(), sa.as_mut_ptr(), 1, 0) },
-1
);
assert_eq!(
unsafe { divsufsort(text.as_ptr(), ptr::null_mut(), 1, 0) },
-1
);
assert_eq!(
unsafe { divsufsort(text.as_ptr(), sa.as_mut_ptr(), -1, 0) },
-1
);
}
#[test]
fn sorts_trivial_inputs() {
/* empty */
let text = [0u8; 1];
let mut sa = [i32::MIN; 1];
assert_eq!(
unsafe { divsufsort(text.as_ptr(), sa.as_mut_ptr(), 0, 0) },
0
);
assert_eq!(sa[0], i32::MIN, "n == 0 must not touch the output");
/* single byte */
assert_eq!(build_sa(b"z"), [0]);
/* two bytes: ascending, descending, and equal */
assert_eq!(build_sa(b"ab"), [0, 1]);
assert_eq!(build_sa(b"ba"), [1, 0]);
assert_eq!(build_sa(b"aa"), [1, 0]);
}
#[test]
fn sorts_all_equal_bytes() {
let text = vec![b'q'; 10_000];
let sa = build_sa(&text);
/* For a constant text the shortest suffix sorts first. */
let expected: Vec<i32> = (0..text.len() as i32).rev().collect();
assert_eq!(sa, expected);
}
#[test]
fn sorts_abracadabra_exactly() {
/* Hand-computed: a(10) abra(7) abracadabra(0) acadabra(3) adabra(5)
* bra(8) bracadabra(1) cadabra(4) dabra(6) ra(9) racadabra(2). */
assert_eq!(build_sa(b"abracadabra"), [10, 7, 0, 3, 5, 8, 1, 4, 6, 9, 2]);
}
#[test]
fn matches_reference_on_periodic_text() {
/* Tandem repeats exercise trsort's repeat partitioning. */
let text: Vec<u8> = b"ab".iter().copied().cycle().take(4096).collect();
let sa = build_sa(&text);
assert_valid_sa(&text, &sa);
assert_eq!(sa, reference_sa(&text));
}
#[test]
fn matches_reference_on_random_bytes() {
let text = lcg_bytes(8192, 0x0BAD_5EED, 256);
let sa = build_sa(&text);
assert_valid_sa(&text, &sa);
assert_eq!(sa, reference_sa(&text));
}
#[test]
fn matches_reference_on_low_alphabet_text() {
/* A four-symbol alphabet produces the large first-two-character
* buckets that reach sssort's block merging and the deeper trsort
* paths. */
let text = lcg_bytes(16_384, 0xDEAD_BEEF, 4);
let sa = build_sa(&text);
assert_valid_sa(&text, &sa);
assert_eq!(sa, reference_sa(&text));
}
}
+2
View File
@@ -5,6 +5,8 @@ pub mod bitstream;
pub mod common; pub mod common;
pub mod cpu; pub mod cpu;
pub mod debug; pub mod debug;
#[cfg(feature = "dict-builder")]
pub mod divsufsort;
pub mod entropy_common; pub mod entropy_common;
pub mod errors; pub mod errors;
#[cfg(feature = "compression")] #[cfg(feature = "compression")]
+2
View File
@@ -79,6 +79,8 @@ RUST_TARGET_DIR := $(RUST_DIR)/target/$(RUST_BUILD_CONFIG)
RUST_STATICLIB := $(RUST_TARGET_DIR)/release/libzstd_rs.a RUST_STATICLIB := $(RUST_TARGET_DIR)/release/libzstd_rs.a
RUST_TARGET_32 ?= i686-unknown-linux-gnu RUST_TARGET_32 ?= i686-unknown-linux-gnu
RUST_STATICLIB_32 := $(RUST_TARGET_DIR)/$(RUST_TARGET_32)/release/libzstd_rs.a RUST_STATICLIB_32 := $(RUST_TARGET_DIR)/$(RUST_TARGET_32)/release/libzstd_rs.a
# Tests build every library module, so they use the crate's default feature
# set (compression, decompression, and dict-builder) plus any forced HUF mode.
RUST_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \ RUST_CARGO_FLAGS := --manifest-path $(RUST_MANIFEST) --release \
--target-dir $(RUST_TARGET_DIR) --target-dir $(RUST_TARGET_DIR)
ifneq ($(RUST_HUF_FEATURE),) ifneq ($(RUST_HUF_FEATURE),)