]> Cypherpunks repositories - gostls13.git/commitdiff
suffixarray: faster creation algorithm
authorEric Eisner <eric.d.eisner@gmail.com>
Wed, 12 Jan 2011 05:46:50 +0000 (21:46 -0800)
committerRobert Griesemer <gri@golang.org>
Wed, 12 Jan 2011 05:46:50 +0000 (21:46 -0800)
This implements the algorithm qsufsort using the sort package
as a sorting primitive. Its worst-case performance is O(N*log(N)), and it
uses only an additional slice of N ints of memory during creation.

Benchmarks (seconds):
           old    new
10k nulls          149    0.044
1M English corpus  32.0   3.6

R=gri, gri1
CC=golang-dev
https://golang.org/cl/3752044

src/pkg/index/suffixarray/Makefile
src/pkg/index/suffixarray/qsufsort.go [new file with mode: 0644]
src/pkg/index/suffixarray/suffixarray.go
src/pkg/index/suffixarray/suffixarray_test.go

index 626ec406aee52bead078e59c5190309316054438..297c4279f929192cc53cb43f31fbeef72ee8157b 100644 (file)
@@ -6,6 +6,7 @@ include ../../../Make.inc
 
 TARG=index/suffixarray
 GOFILES=\
+       qsufsort.go\
        suffixarray.go\
 
 include ../../../Make.pkg
diff --git a/src/pkg/index/suffixarray/qsufsort.go b/src/pkg/index/suffixarray/qsufsort.go
new file mode 100644 (file)
index 0000000..0e6894a
--- /dev/null
@@ -0,0 +1,164 @@
+// Copyright 2011 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+// This algorithm is based on "Faster Suffix Sorting"
+//   by N. Jesper Larsson and Kunihiko Sadakane
+// paper: http://www.larsson.dogma.net/ssrev-tr.pdf
+// code:  http://www.larsson.dogma.net/qsufsort.c
+
+// This algorithm computes the suffix array sa by computing its inverse.
+// Consecutive groups of suffixes in sa are labeled as sorted groups or
+// unsorted groups. For a given pass of the sorter, all suffixes are ordered
+// up to their first h characters, and sa is h-ordered. Suffixes in their
+// final positions and unambiguouly sorted in h-order are in a sorted group.
+// Consecutive groups of suffixes with identical first h characters are an
+// unsorted group. In each pass of the algorithm, unsorted groups are sorted
+// according to the group number of their following suffix.
+
+// In the implementation, if sa[i] is negative, it indicates that i is
+// the first element of a sorted group of length -sa[i], and can be skipped.
+// An unsorted group sa[i:k] is given the group number of the index of its
+// last element, k-1. The group numbers are stored in the inverse slice (inv),
+// and when all groups are sorted, this slice is the inverse suffix array.
+
+package suffixarray
+
+import "sort"
+
+func qsufsort(data []byte) []int {
+       // initial sorting by first byte of suffix
+       sa := sortedByFirstByte(data)
+       if len(sa) < 2 {
+               return sa
+       }
+       // initialize the group lookup table
+       // this becomes the inverse of the suffix array when all groups are sorted
+       inv := initGroups(sa, data)
+
+       // the index starts 1-ordered
+       sufSortable := &suffixSortable{sa, inv, 1}
+
+       for sa[0] > -len(sa) { // until all suffixes are one big sorted group
+               // The suffixes are h-ordered, make them 2*h-ordered
+               pi := 0 // pi is first position of first group
+               sl := 0 // sl is negated length of sorted groups
+               for pi < len(sa) {
+                       if s := sa[pi]; s < 0 { // if pi starts sorted group
+                               pi -= s // skip over sorted group
+                               sl += s // add negated length to sl
+                       } else { // if pi starts unsorted group
+                               if sl != 0 {
+                                       sa[pi+sl] = sl // combine sorted groups before pi
+                                       sl = 0
+                               }
+                               pk := inv[s] + 1 // pk-1 is last position of unsorted group
+                               sufSortable.sa = sa[pi:pk]
+                               sort.Sort(sufSortable)
+                               sufSortable.updateGroups(pi)
+                               pi = pk // next group
+                       }
+               }
+               if sl != 0 { // if the array ends with a sorted group
+                       sa[pi+sl] = sl // combine sorted groups at end of sa
+               }
+
+               sufSortable.h *= 2 // double sorted depth
+       }
+
+       for i := range sa { // reconstruct suffix array from inverse
+               sa[inv[i]] = i
+       }
+       return sa
+}
+
+
+func sortedByFirstByte(data []byte) []int {
+       // total byte counts
+       var count [256]int
+       for _, b := range data {
+               count[b]++
+       }
+       // make count[b] equal index of first occurence of b in sorted array
+       sum := 0
+       for b := range count {
+               count[b], sum = sum, count[b]+sum
+       }
+       // iterate through bytes, placing index into the correct spot in sa
+       sa := make([]int, len(data))
+       for i, b := range data {
+               sa[count[b]] = i
+               count[b]++
+       }
+       return sa
+}
+
+
+func initGroups(sa []int, data []byte) []int {
+       // label contiguous same-letter groups with the same group number
+       inv := make([]int, len(data))
+       prevGroup := len(sa) - 1
+       groupByte := data[sa[prevGroup]]
+       for i := len(sa) - 1; i >= 0; i-- {
+               if b := data[sa[i]]; b < groupByte {
+                       if prevGroup == i+1 {
+                               sa[i+1] = -1
+                       }
+                       groupByte = b
+                       prevGroup = i
+               }
+               inv[sa[i]] = prevGroup
+               if prevGroup == 0 {
+                       sa[0] = -1
+               }
+       }
+       // Separate out the final suffix to the start of its group.
+       // This is necessary to ensure the suffix "a" is before "aba"
+       // when using a potentially unstable sort.
+       lastByte := data[len(data)-1]
+       s := -1
+       for i := range sa {
+               if sa[i] >= 0 {
+                       if data[sa[i]] == lastByte && s == -1 {
+                               s = i
+                       }
+                       if sa[i] == len(sa)-1 {
+                               sa[i], sa[s] = sa[s], sa[i]
+                               inv[sa[s]] = s
+                               sa[s] = -1 // mark it as an isolated sorted group
+                               break
+                       }
+               }
+       }
+       return inv
+}
+
+
+type suffixSortable struct {
+       sa  []int
+       inv []int
+       h   int
+}
+
+func (x *suffixSortable) Len() int           { return len(x.sa) }
+func (x *suffixSortable) Less(i, j int) bool { return x.inv[x.sa[i]+x.h] < x.inv[x.sa[j]+x.h] }
+func (x *suffixSortable) Swap(i, j int)      { x.sa[i], x.sa[j] = x.sa[j], x.sa[i] }
+
+
+func (x *suffixSortable) updateGroups(offset int) {
+       prev := len(x.sa) - 1
+       group := x.inv[x.sa[prev]+x.h]
+       for i := prev; i >= 0; i-- {
+               if g := x.inv[x.sa[i]+x.h]; g < group {
+                       if prev == i+1 { // previous group had size 1 and is thus sorted
+                               x.sa[i+1] = -1
+                       }
+                       group = g
+                       prev = i
+               }
+               x.inv[x.sa[i]] = prev + offset
+               if prev == 0 { // first group has size 1 and is thus sorted
+                       x.sa[0] = -1
+               }
+       }
+}
index 88cf925fc2bcabda72e75d55f1d6a5a00cc9407f..628e000e1d35a5a3601f09ad0db7c2ae4b71b0e0 100644 (file)
@@ -22,11 +22,6 @@ import (
        "sort"
 )
 
-// BUG(gri): For larger data (10MB) which contains very long (say 100000)
-// contiguous sequences of identical bytes, index creation time will be extremely slow.
-
-// TODO(gri): Use a more sophisticated algorithm to create the suffix array.
-
 
 // Index implements a suffix array for fast substring search.
 type Index struct {
@@ -36,16 +31,9 @@ type Index struct {
 
 
 // New creates a new Index for data.
-// Index creation time is approximately O(N*log(N)) for N = len(data).
-//
+// Index creation time is O(N*log(N)) for N = len(data).
 func New(data []byte) *Index {
-       sa := make([]int, len(data))
-       for i := range sa {
-               sa[i] = i
-       }
-       x := &Index{data, sa}
-       sort.Sort((*index)(x))
-       return x
+       return &Index{data, qsufsort(data)}
 }
 
 
@@ -192,12 +180,3 @@ func (x *Index) FindAllIndex(r *regexp.Regexp, n int) (result [][]int) {
        }
        return
 }
-
-
-// index is used to hide the sort.Interface
-type index Index
-
-func (x *index) Len() int           { return len(x.sa) }
-func (x *index) Less(i, j int) bool { return bytes.Compare(x.at(i), x.at(j)) < 0 }
-func (x *index) Swap(i, j int)      { x.sa[i], x.sa[j] = x.sa[j], x.sa[i] }
-func (a *index) at(i int) []byte    { return a.data[a.sa[i]:] }
index 659bce042c084261768536fe938e6f6d6f5972fb..b3486a96d08a9fabf3722ab770610159c697c498 100644 (file)
@@ -5,6 +5,7 @@
 package suffixarray
 
 import (
+       "bytes"
        "container/vector"
        "regexp"
        "sort"
@@ -204,9 +205,26 @@ func testLookups(t *testing.T, tc *testCase, x *Index, n int) {
 }
 
 
+// index is used to hide the sort.Interface
+type index Index
+
+func (x *index) Len() int           { return len(x.sa) }
+func (x *index) Less(i, j int) bool { return bytes.Compare(x.at(i), x.at(j)) < 0 }
+func (x *index) Swap(i, j int)      { x.sa[i], x.sa[j] = x.sa[j], x.sa[i] }
+func (a *index) at(i int) []byte    { return a.data[a.sa[i]:] }
+
+
+func testConstruction(t *testing.T, tc *testCase, x *Index) {
+       if !sort.IsSorted((*index)(x)) {
+               t.Errorf("testConstruction failed %s", tc.name)
+       }
+}
+
+
 func TestIndex(t *testing.T) {
        for _, tc := range testCases {
                x := New([]byte(tc.source))
+               testConstruction(t, &tc, x)
                testLookups(t, &tc, x, 0)
                testLookups(t, &tc, x, 1)
                testLookups(t, &tc, x, 10)