This chapter introduces the linear-time in-place algorithm for suffix sorting1, Optimal In-Place Suffix Sorting.
Global setup
The target string is , and the suffix array is . String indices start at 0, and the final character is a sentinel, which we may assume is 0.
Suffix sorting over an integer alphabet
This part can in fact be viewed as an in-place version of SA-IS.
It is the algorithm whose details are clearest in the original paper and whose implementation is relatively simple. It also forms the basis for understanding the later algorithms, so it is the focus of this article.
The in-place principle is to use a renamed in place of the S and L buckets, and to use additional operations in place of the type buckets.
Renaming the target string Pat
In brief, rename without changing the relative order of the suffixes, then use the renamed in place of the original S and L buckets to indicate bucket heads and tails.
Rename each S-type character in to the index of its bucket tail, and each L-type character to the index of its bucket head.
The process is illustrated below:
The renamed (from now on, call the renamed simply ) is:
Within a bucket, L-type suffixes are smaller and serve as the bucket head, while S-type suffixes are larger and serve as the bucket tail, so the relative suffix order is preserved.
The renaming steps are:
- As in SA-IS, count every character in and compute prefix sums (counting sort) to build the S/L buckets; here, however, store the prefix sums in .
- Scan every character of from right to left. By recording only the type of the previous character, determine each character's type online and rename it according to the prefix sums.
Sorting LMS characters
The key idea here is an internal-counter technique.
Initialization
Initially set every entry of to E (EMPTY).
Scan from right to left. Whenever an LMS character is found, set the marker at :
If is E, set it to U (UNIQUE).
If is U, set it to M (MULTIPLE).
Do nothing in all other cases.
The result is shown below:
Placing LMS-character indices in SA
Scan from right to left. For an LMS character , distinguish cases according to the marker in :
U: directly set .
M: this means that the bucket contains at least two LMS characters.
-
If the bucket contains at least three LMS characters:
- Use the penultimate bucket position as a temporary counter recording how many LMS characters in the bucket have been filled (the last position carries marker M).
- Insert each new LMS character starting at the third position from the end, incrementing the temporary counter.
- When the bucket becomes full, shift every element from the bucket head through the third-from-last position two places to the right, then insert the new element at the second bucket position (set the first position to E).
-
If the bucket contains exactly two LMS characters, no counter is needed; insert them directly from right to left.
Normal values:
According to the discussion above, whether a bucket contains two or more LMS characters, this means that is the position of the last LMS character still to be inserted into the bucket.
Scan left from the bucket head until finding the first position marked E, and set it to .
Finally, scan from right to left to clear any remaining special marker M (the bucket was not full, so M and the counter were not overwritten).
Shift the bucket's LMS characters two positions to the right as above, and fill the positions vacated on the left with E.
The result is:
Because each bucket is moved and scanned only once in this phase, the time complexity is .
Induced sorting of LMS substrings
Induced sorting of LMS prefixes
Induced-sorting LMS prefixes is exactly the same as induced-sorting suffixes later in the SA-IS process (using the same function), so we omit it here.
The sorted result is:
Moving sorted LMS substrings to the end of SA
Building the reduced subproblem string Pat1
Scan the LMS substrings at the end of from left to right, determine their relative order (“rename”), and store the renamed value of in .
Because LMS characters are not adjacent, there is no collision. This places the renamed values in the original order of the substrings they represent in :
Then scan and collect these renamed values at the head of :
Solving Pat1 recursively to finish sorting LMS suffixes
As in SA-IS, recursively sort the suffixes of the reduced at the head of , and store the result at the end of :
Move the tail's to the head of , then scan from right to left again, place its LMS characters in their order in at the end of :
Using the “lookup table” at the end of , restore the head values of to the index positions of the corresponding LMS suffixes in :
Place the sorted LMS suffixes from the head of into their corresponding buckets in order, starting from the bucket tails:
Induced sorting of all suffixes in Pat1
This part uses the internal-counter technique introduced earlier to perform in-place induced sorting.
Suppose the sorted LMS suffixes are already at the bucket tails. We use them to induce L-type suffixes4:
As when sorting LMS characters, first count L-type characters with special markers:
Scan SA from left to right. As with LMS sorting, the more complicated part is determining the type of , which requires case analysis (see the code for details):
Unlike SA-IS, after inducing one character type, we must clear the LMS characters to avoid interfering with the next in-place induced sort:
Inducing S-type suffixes from L-type suffixes is completely symmetric to inducing L-type suffixes from LMS suffixes, so we do not repeat the details.
At this point, induced sorting is complete.
Implementation
In terms of running time, this has no significant difference from SA-IS. Its space usage is less than of the original (at the cost of twice the code), which makes the title Optimal In-Place Suffix Sorting1 well deserved.
Reference code
use std::cmp::max;
use std::cmp::Ordering;
use std::slice::from_raw_parts_mut;
const LTYPE: bool = false;
const STYPE: bool = true;
const MAX_SA_VALUE: usize = usize::MAX / 2;
const EMPTY: usize = MAX_SA_VALUE + 1;
const UNIQUE: usize = MAX_SA_VALUE + 2;
const MULTI: usize = MAX_SA_VALUE + 3; // >= 258
fn lms_str_cmp<E: Ord>(l1: &[E], l2: &[E]) -> Ordering {
for (x, y) in l1.iter().zip(l2.iter()) {
let cmp_res = x.cmp(&y);
if cmp_res != Ordering::Equal { return cmp_res; }
}
Ordering::Equal
}
#[inline]
fn pat_char_type(cur: usize, prev: usize, last_scanned_type: bool) -> bool {
if cur < prev || cur == prev && last_scanned_type == STYPE { STYPE }
else { LTYPE }
}
fn rename_pat(pat: &mut [usize], sa: &mut [usize]) {
let patlastpos = pat.len() - 1;
// Convert everything to bucket heads
//sa.fill(0);
for i in 0..sa.len() { sa[i] = 0 }
for i in 0..pat.len() { sa[pat[i]] += 1 }
for i in 1..sa.len() { sa[i] += sa[i - 1] }
for i in 0..pat.len() - 1 {
pat[i] = sa[pat[i]] - 1;
};
// Convert L-suffixes to bucket heads
//sa.fill(0);
for i in 0..sa.len() { sa[i] = 0 }
for i in 0..pat.len() { sa[pat[i]] += 1 }
let mut last_scanned_type = STYPE;
pat[patlastpos] = 0;
for i in (0..pat.len() - 1).rev() {
if pat_char_type(pat[i], pat[i + 1], last_scanned_type) == STYPE {
last_scanned_type = STYPE;
} else {
pat[i] -= sa[pat[i]] - 1;
last_scanned_type = LTYPE;
}
}
}
fn sort_lms_char(pat: &mut [usize], sa: &mut [usize]) -> usize {
//sa.fill(EMPTY);
for i in 0..sa.len() { sa[i] = EMPTY }
let mut last_scanned_type = STYPE;
for i in (0..pat.len() - 1).rev() {
if pat_char_type(pat[i], pat[i + 1], last_scanned_type) == STYPE {
last_scanned_type = STYPE;
} else {
if last_scanned_type == STYPE { // pat[i + 1] is LMS type
sa[pat[i + 1]] += 1;
}
last_scanned_type = LTYPE;
}
}
let mut lms_cnt = 0;
last_scanned_type = STYPE;
for i in (0..pat.len() - 1).rev() {
if pat_char_type(pat[i], pat[i + 1], last_scanned_type) == STYPE {
last_scanned_type = STYPE;
} else {
let e_i = i + 1;
let e = pat[e_i];
if last_scanned_type == STYPE { // pat[i + 1] is LMS type
lms_cnt += 1;
if sa[e] == UNIQUE {
sa[e] = e_i;
} else if sa[e] >= MULTI && sa[e - 1] == EMPTY {
if sa[e - 2] == EMPTY {
sa[e - 2] = e_i;
sa[e - 1] = 1; // set counter
} else { // MUL = 2
sa[e] = e_i;
sa[e - 1] = EMPTY;
}
} else if sa[e] >= MULTI && sa[e - 1] != EMPTY {
let c = sa[e - 1]; // get counter
if sa[e - 2 - c] == EMPTY {
sa[e - 2 - c] = e_i;
sa[e - 1] += 1; // update counter
} else {
for j in (1..c + 1).rev() {
sa[e - c + j] = sa[e - 2 - c + j]
}
sa[e - c] = e_i;
sa[e - c - 1] = EMPTY;
}
} else if sa[e] < EMPTY {
for j in (0..e).rev() {
if sa[j] == EMPTY {
sa[j] = e_i;
break;
}
}
}
}
last_scanned_type = LTYPE;
}
}
for i in (0..pat.len()).rev() {
if sa[i] >= MULTI {
let c = sa[i - 1];
for j in (1..c + 1).rev() { // Reverse order prevents earlier entries from overwriting later ones
sa[i - c + j] = sa[i - 2 - c + j];
}
sa[i - c - 1] = EMPTY;
sa[i - c] = EMPTY;
}
}
lms_cnt
}
fn sort_lms_substr(pat: &mut [usize], sa: &mut [usize]) {
// step 1
induced_sort(pat, sa);
// step 2
let pat_last_pos = pat.len() - 1;
let mut lms_cnt = 0;
let mut i = pat_last_pos;
let mut bucket_tail_ptr = pat_last_pos + 1; // for renamed bucket ver
let mut bucket = EMPTY; // This can be omitted, but it makes the code easier to write
let mut num = 0; // S type number of bucket
while i > 0 {
if pat[sa[i]] != bucket { // reach new bucket
num = 0;
let mut l = 0;
while pat[sa[i - l]] == pat[sa[i]] { // Scan the bucket to count S characters; by definition, the loop ends when l=i
let pat_i = sa[i - l]; // l < i, so i - l > 0 and 0 <= pat_i < patlen - 1
if pat[pat_i] < pat[pat_i + 1] {
let mut k = pat_i;
while k > 0 && pat[k - 1] == pat[pat_i] { k -= 1 }
num += pat_i - k + 1;
} else {
break; // The bucket contains no S characters; stop scanning
}
l += 1;
}
bucket_tail_ptr = i;
bucket = pat[sa[bucket_tail_ptr]];
}
if num > 0
&& i > bucket_tail_ptr - num
&& sa[i] > 0
&& pat[sa[i]] < pat[sa[i] - 1] {
sa[pat_last_pos - lms_cnt] = sa[i];
lms_cnt += 1;
}
i -= 1;
}
sa[pat_last_pos - lms_cnt ] = sa[i]; // i = 0
lms_cnt += 1;
//sa[0..pat_last_pos - lms_cnt + 1].fill(EMPTY);
for i in 0..pat_last_pos - lms_cnt + 1 { sa[i] = EMPTY }
}
fn construct_pat1(pat: &mut [usize], sa: &mut [usize], lms_cnt: usize) -> bool {
let patlen = pat.len();
let mut prev_lms_str_len = 1;
let mut rank = 0;
sa[(patlen - 1) / 2] = rank;
let mut has_duplicated_char = false;
for i in patlen - lms_cnt + 1..patlen { // Start from the character after the sentinel
let mut j = sa[i];
while pat[j] <= pat[j + 1] { j += 1 } // Find the first L character to the right of suf(sa[i]); excluding the sentinel LMS suffix prevents overflow
let mut k = j;
while k + 1 < patlen && pat[k] >= pat[k + 1] { k += 1 } // Find the first LMS character to the right of suf(sa[i])
let cur_lms_str_len = k + 1 - sa[i];
let cmp_res = lms_str_cmp(&pat[sa[i]..sa[i] + cur_lms_str_len], &pat[sa[i - 1]..sa[i - 1] + prev_lms_str_len]);
if cmp_res != Ordering::Equal {
rank += 1
}
if rank == sa[sa[i - 1] / 2] {
has_duplicated_char = true;
}
let rank_index = sa[i] / 2;
sa[rank_index] = rank; // Integer division
prev_lms_str_len = cur_lms_str_len;
}
// move to head of sa
let mut j = 0;
for i in 0..patlen - lms_cnt {
if sa[i] != EMPTY {
sa[j] = sa[i];
if i > j {
sa[i] = EMPTY;
}
j += 1;
}
}
//sa[lms_cnt..patlen].fill(EMPTY);
for i in lms_cnt..patlen { sa[i] = EMPTY }
has_duplicated_char
}
fn sort_lms_suf(pat: &mut [usize], sa: &mut [usize], lms_cnt: usize, has_duplicated_char: bool) {
// solve T1 recursively
let patlen = pat.len();
let salen = sa.len();
unsafe {
let sa_ptr = sa.as_mut_ptr();
let mut pat1 = from_raw_parts_mut(sa_ptr, lms_cnt);
let mut sa1 = from_raw_parts_mut(sa_ptr.offset((patlen - lms_cnt) as isize), salen - (patlen - lms_cnt));
if has_duplicated_char {
_compute_suffix_array_16_1(&mut pat1, &mut sa1);
} else {
for i in 0..lms_cnt { sa1[pat1[i]] = i }
}
}
// move SA1 to SA[0...n1-1]
for i in 0..lms_cnt {
sa[i] = sa[patlen- lms_cnt + i];
}
// put all LMS-suffixes in SA tail
let mut last_scanned_type = STYPE;
let mut j = 0;
for i in (0..pat.len() - 1).rev() {
if pat[i] < pat[i + 1] || pat[i] == pat[i + 1] && last_scanned_type == STYPE {
last_scanned_type = STYPE;
} else {
if last_scanned_type == STYPE {
sa[patlen - 1 - j] = i + 1;
j += 1;
}
last_scanned_type = LTYPE;
}
}
// backward map the LMS-suffixes rank
for i in 0..lms_cnt {
let relative_rank = sa[i];
sa[i] = sa[patlen - lms_cnt + relative_rank];
sa[patlen - lms_cnt + relative_rank] = EMPTY;
}
let mut tail = EMPTY;
let mut rfp = EMPTY;
for i in (1..lms_cnt).rev() { // Keep sa[0] in place
if pat[sa[i]] != tail {
tail = pat[sa[i]];
rfp = tail;
}
sa[rfp] = sa[i];
if rfp != i { sa[i] = EMPTY }
rfp -= 1;
}
}
// PASS!
fn induced_sort(pat: &mut [usize], sa: &mut [usize]) {
let patlen = pat.len();
// place L-suff in SA
// init
let mut last_scanned_type = STYPE;
for i in (0..patlen - 1).rev() {
if pat_char_type(pat[i], pat[i + 1], last_scanned_type) == LTYPE {
sa[pat[i]] += 1; // >= EMPTY
last_scanned_type = LTYPE;
} else {
last_scanned_type = STYPE;
}
}
//place
let mut i = 0;
while i < patlen {
if sa[i] < EMPTY && sa[i] > 0 {
let j = sa[i] - 1;
let mut is_ltype = false;
if pat[j] > pat[j + 1] {
is_ltype = true;
} else if pat[j] == pat[j + 1] { // Determine whether sa[i] is the index of an L-suffix
let next_i = sa[pat[sa[i]]];
if next_i >= MULTI {
is_ltype = true;
} else if next_i < EMPTY && pat[sa[i]] + 1 < patlen {
if sa[pat[sa[i]] + 1] == EMPTY {
is_ltype = true;
} else if sa[pat[sa[i]] + 1] < EMPTY {
if pat[sa[pat[sa[i]] + 1]] == pat[sa[i]] {
is_ltype = true;
}
}
}
}
if is_ltype {
if sa[pat[j]] == UNIQUE {
sa[pat[j]] = j;
} else if sa[pat[j]] >= MULTI && sa[pat[j] + 1] == EMPTY {
if sa[pat[j]] - EMPTY > 2 {
sa[pat[j] + 2] = j;
sa[pat[j] + 1] = 1; // set counter
} else {
sa[pat[j]] = j;
}
} else if sa[pat[j]] >= MULTI && sa[pat[j] + 1] != EMPTY {
let e = pat[j];
let c = sa[e + 1];
let lfp = e + c + 2;
if c + 2 < sa[pat[j]] - EMPTY { // The bucket tail has not been reached
sa[lfp] = j;
sa[e + 1] += 1; // update counter
} else {
for k in 1..c + 1 {
sa[e + k - 1] = sa[e + k + 1];
}
sa[e + c] = j;
sa[e + c + 1] = EMPTY;
if i >= e + 2 && i <= e + c + 1 {
i -= 2;
}
}
} else if sa[pat[j]] < EMPTY {
for k in pat[j]..patlen {
if sa[k] == EMPTY {
sa[k] = j;
break;
}
}
}
}
} else if sa[i] >= MULTI {
i += 1;
}
i += 1;
}
// Remove LMS suffixes from SA; a bucket may contain multiple LMS suffixes
last_scanned_type = STYPE;
for i in (0..pat.len() - 1).rev() {
if pat_char_type(pat[i], pat[i + 1], last_scanned_type) == STYPE {
last_scanned_type = STYPE;
} else {
if last_scanned_type == STYPE { // pat[i + 1] is LMS type
if sa[pat[i + 1]] <= EMPTY {
sa[pat[i + 1]] = UNIQUE;
} else {
sa[pat[i + 1]] += 1;
}
}
last_scanned_type = LTYPE;
}
}
i = patlen - 1;
while i > 0 {
if sa[i] > EMPTY {
let c = sa[i] - EMPTY;
for k in 0..c {
sa[i - k] = EMPTY;
}
i -= c - 1;
}
i -= 1;
}
sa[0] = pat.len() - 1;
// place S-suff in SA
// init
let mut last_scanned_type = STYPE;
for i in (0..patlen - 1).rev() {
if pat_char_type(pat[i], pat[i + 1], last_scanned_type) == STYPE {
if sa[pat[i]] >= EMPTY {
sa[pat[i]] += 1;
} else {
sa[pat[i]] = UNIQUE;
}
last_scanned_type = STYPE;
} else {
last_scanned_type = LTYPE;
}
}
i = patlen - 1;
while i > 0 {
if sa[i] < EMPTY && sa[i] > 0 {
let j = sa[i] - 1;
let mut is_stype = false;
if pat[j] < pat[j + 1] {
is_stype = true;
} else if pat[j] == pat[j + 1] { // Determine whether sa[i] is the index of an S-suffix
let next_i = sa[pat[sa[i]]];
if next_i >= MULTI {
is_stype = true;
} else if next_i < EMPTY && pat[sa[i]] - 1 > 0 {
if sa[pat[sa[i]] - 1] == EMPTY {
is_stype = true;
} else if sa[pat[sa[i]] - 1] < EMPTY {
if pat[sa[pat[sa[i]] - 1]] == pat[sa[i]] {
is_stype = true;
}
}
}
}
if is_stype {
if sa[pat[j]] == UNIQUE {
sa[pat[j]] = j;
} else if sa[pat[j]] >= MULTI && sa[pat[j] - 1] == EMPTY {
if sa[pat[j]] - EMPTY > 2 {
sa[pat[j] - 2] = j;
sa[pat[j] - 1] = 1; // set counter
} else {
sa[pat[j]] = j;
}
} else if sa[pat[j]] >= MULTI && sa[pat[j] - 1] != EMPTY {
let e = pat[j];
let c = sa[e - 1];
let num = sa[pat[j]] - EMPTY;
if c + 2 < num { // The bucket head has not been reached
let rfp = e - c - 2;
sa[rfp] = j;
sa[e - 1] += 1;
} else {
for k in 1..c + 1 {
sa[e - k + 1] = sa[e - k - 1];
}
sa[e - c] = j;
sa[e - c - 1] = EMPTY;
if i >= e - num + 1 && i <= e - 2 {
i += 2;
}
}
} else if sa[pat[j]] < EMPTY {
for k in (0..pat[j]).rev() {
if sa[k] == EMPTY {
sa[k] = j;
break;
}
}
}
}
} else if sa[i] >= MULTI {
i -= 1;
}
i -= 1;
}
}
fn _compute_suffix_array_16_1(pat: &mut [usize], sa: &mut [usize]) {
rename_pat(pat, sa);
let lms_cnt = sort_lms_char(pat, sa);
sort_lms_substr(pat, sa);
let has_duplicated_char = construct_pat1(pat, sa, lms_cnt);
sort_lms_suf(pat, sa, lms_cnt, has_duplicated_char);
induced_sort(pat, sa);
}
pub fn suffix_array_16(pat: &[u8]) -> Vec<usize> {
let mut pat = pat.into_iter().map(|x| *x as usize).collect::<Vec<usize>>();
pat.push(0);
let mut sa = vec![0; max(pat.len(), 256) * 1];
_compute_suffix_array_16_1(&mut pat[..], &mut sa[..]);
sa
}
fn input() -> String {
use std::io;
let mut input = String::new();
io::stdin().read_line(&mut input).unwrap();
String::from(input.trim())
}
fn main() {
let pat = input();
let sa_16 = suffix_array_16(pat.as_bytes());
for i in 1..pat.len() + 1 { print!("{} ", sa_16[i] + 1) }
}
Suffix sorting over an integer alphabet in read-only memory
Use a sophisticated method to solve a sophisticated problem: use divide and conquer to handle tight space limits.
The challenge in the implementation is to build BitMaps on 5, replacing the bucket-head and bucket-tail positions originally indicated by the renamed T.
Here, BitMaps means an ordered dictionary (multiset) represented by a bit vector. It is a compact data structure.
For now, anyone interested must consult the original paper and the papers on BitMaps cited by this article.
Suffix sorting over a general alphabet in read-only memory
The prerequisites are merge sort and heap sort.
Because the author has doubts about the time complexity of determining character types in this method, it is not described here; consult the original paper for details.
Notes
Footnotes
-
Li, Zhize; Li, Jian; Huo, Hongwei (2016). Optimal In-Place Suffix Sorting. Proceedings of the 25th International Symposium on String Processing and Information Retrieval (SPIRE). Lecture Notes in Computer Science. 1114. ↩ ↩2
-
Ge Nong, Sen Zhang, and Wai Hong Chan. Linear suffix array construction by almost pure induced-sorting. In Data Compression Conference (DCC), pages 193–202. IEEE, 2009. ↩
-
We recommend reading this blog post and its issue list. ↩
-
For an LML suffix, induce S-type suffixes first. The only difference is that the sentinel must also be included when computing an LML suffix. ↩
-
Gonzalo Navarro and Eliana Providel. Fast, small, simple rank/select on bitmaps. In Proceedings of the 11th International Symposium on Experimental Algorithms (SEA), pages 295–306, 2012. ↩