SA-IS 学习笔记

发布时间 2023-10-02 20:32:03作者: 383494

目录

  1. 鲜花

  2. 约定

  3. 后缀类型

  4. 诱导排序

  5. 算法过程

  6. 性能测试

  7. 推荐阅读

闲话

我太蒻了,学了两三天才会。/kel

好像也没什么 duliu 出题人卡 \(O(n \log n)\) 的 SA...?

网上大多的 blog 都说要在字符串尾部加一个字典序为 \(+\infty\) 的字符,但这样感觉不符合部分定义且不符合习惯,因此我加的是 \(-\infty\)\0.

本文中大多数定理都不证,因为读者自证不难我不会,且可以在其他 blog 中找到证明。

可能后面我学会了之后会来补证明.

约定

  • 字符串下标从 \(0\) 开始。

  • 用大写字母表示字符串,小写字母表示字符/数字。

  • 对于一个字符串 \(S\) 的后缀 \(A\)\(B\),可以用它们的起始位置 \(a\)\(b\)\(S^{a}, S^{b}\) 代指,一般从上下文中可以看出。不用 \(S_a,S_b\) 的原因是它们用来代表单个字符。

后缀类型

对于 \(a < a+1, S^{a} < S^{a+1}\) 的串 \(A\),称其为 Small(S)型。

否则,称其为 Large(L)型。

特别地,对于左边为 L 型的 S 型,又可称为 Left-Most-Small(LMS)型。

记后缀 \(A\) 的类型为 \(\mathrm{typ}(A) = \mathrm{typ}(a) \in \{S, L\}\)

定理:\(S_a = S_b \implies \mathrm{typ}(a) = \mathrm{typ}(b)\).

从后向前扫一遍可以在 \(\Theta(n)\) 的时间内得出 \(S\)\(\mathrm{typ}\) 数组。

诱导排序

  1. 预处理出每种首字母的出现次数,不难发现对于某种特定的首字母,以它开头的后缀一定落在 SA 数组某个区间内。这个区间可以通过计算(出现次数)的前缀和得到。

  2. 将 SA 数组填满 \(-1\).

  3. 将所有 LMS 型后缀按某种顺序插入 SA 数组中。具体地,若当前插入的是 \(A\),则将其放入 \(A\) 开头的区间末尾,除非末尾已经有另一个 LMS 后缀了,此时向前放一格直到找到空位为止。

  4. 按顺序扫描 SA 数组,若 \(\mathrm{SA}_{i-1} \not = -1 \land \mathrm{typ}(\mathrm{SA}_{i-1}) = L\),则将其放入 \(S_{i-1}\) 开头的区间开头。

  5. 将之前加入 SA 的 LMS 后缀丢弃。可以不动 SA 数组,只将区间末尾的标记复位。

  6. 按倒序扫描 SA 数组,若 \(\mathrm{SA}_{i-1} \not = -1 \land \mathrm{typ}(\mathrm{SA}_{i-1}) = S\),则将其放入 \(S_{i-1}\) 开头的区间末尾。

定理:

若一开始插入 LMS 后缀为字符串中的顺序,则它们在 SA 中的相对顺序为正确的;

若一开始插入 LMS 后缀为 LMS 子串的 SA 数组(注意不是原 \(S\) 的 SA 数组)倒序,则诱导排序后能得到 \(S\) 的 SA 数组。

算法过程

  1. 求出 \(\mathrm{typ}\) 数组

  2. 找出 \(S\) 中的 LMS 开头

  3. 诱导排序

  4. 将 LMS 子串离散化

  5. 若离散后的字符串中每个字符都不同,则基排;否则递归 SA-IS

  6. 诱导排序

代码

#include <iostream>
#include <string>
#include <vector>
#include <algorithm>
#define UP(i,s,e) for(auto i=s; i!=e; ++i)
using std::cin; using std::cout;
constexpr int LEN = 1e6, SIGMA = 128;
namespace SAIS{ // }{{{
// false == L, true == S
template<typename ARRAY>
void sais(ARRAY &str, int len, int sigma, std::vector<int> &sa){ 
    // assert str[len-1] == -INF
    std::vector<bool> typ;
    std::vector<int> bend, lms, bends;
    auto induced_sort = [&](){
        std::vector<int> pos;
        pos.reserve(sigma);
        pos.push_back(0);
        pos.insert(pos.end(), bend.begin(), --bend.end());
        UP(i, 0, len){
            if(sa[i] <= 0) continue;
            if(typ[sa[i]-1] != false) continue;
            sa[pos[str[sa[i]-1]]++] = sa[i]-1;
        }
        pos = bend;
        for(int i=len-1; i>=0; i--){
            if(sa[i] <= 0) continue;
            if(typ[sa[i]-1] != true) continue;
            sa[--pos[str[sa[i]-1]]] = sa[i]-1;
        }
    };
    auto is_lms = [&typ](int x){ if(x == 0) return false; else return typ[x] == true && typ[x-1] == false; };
    auto is_same_lms = [&](int x, int y){
        if(x == -1 || y == -1) return false;
        if(x == y) return true;
        bool moved = false;
        for(;;){
            if(str[x] != str[y]) return false;
            if(typ[x] != typ[y]) return false;
            if(is_lms(x) != is_lms(y)) return false;
            if(is_lms(x) && moved) return true;
            x++; y++;
            moved = true;
        }
    };
    typ.resize(len);
    sa.resize(len);
    bend.resize(sigma);
    std::fill(bend.begin(), bend.end(), 0);
    std::fill(sa.begin(), sa.end(), -1);
    UP(i, 0, len){ bend[str[i]]++; }
    UP(i, 1, sigma){ bend[i] += bend[i-1]; } // lower_bound
    typ[len-1] = true; // S
    for(int i=len-1-1; i>=0; --i){ // get types
        if(str[i] == str[i+1]) typ[i] = typ[i+1];
        else if(str[i] > str[i+1]) typ[i] = false;
        else typ[i] = true;
    }
    UP(i, 1, len){ if(is_lms(i)) lms.push_back(i); }
    bends = bend;
    for(int i:lms){ // insert lms
        sa[--bends[str[i]]] = i;
    }
    induced_sort();
    // discretization
    std::vector<int> lms_names;
    lms_names.resize(len);
    std::fill(lms_names.begin(), lms_names.end(), -1);
    int diffcnt = 0, lastlms = -1; // diffcnt: counter of different lms
    for(int i:sa){
        if(!is_lms(i)) continue;
        if(!is_same_lms(i, lastlms)){
            diffcnt++;
        }
        lastlms = i;
        lms_names[i] = diffcnt-1;
    }
    std::vector<int> nxt_str, nxt_pl;
    nxt_str.reserve(diffcnt);
    nxt_pl.reserve(diffcnt);
    // nxt_sigma == diffcnt
    UP(i, 0, len){
        if(lms_names[i] == -1) continue;
        nxt_str.push_back(lms_names[i]);
        nxt_pl.push_back(i);
    }
    lms_names.clear(); lms_names.shrink_to_fit();
    std::vector<int> sa1;
    if(diffcnt == (int)nxt_str.size()){
        sa1.resize(diffcnt);
        std::fill(sa1.begin(), sa1.end(), -1);
        UP(i, 0, diffcnt){
            sa1[nxt_str[i]] = i;
        }
    } else {
        sais(nxt_str, nxt_str.size(), diffcnt, sa1);
    }
    std::fill(sa.begin(), sa.end(), -1);
    bends = bend;
    for(int i=sa1.size()-1; i>=0; --i){
        int pl = nxt_pl[sa1[i]];
        sa[--bends[str[pl]]] = pl;
    }
    induced_sort();
    return;
}
} // {}}}
namespace m{ // }{{{
std::string str;
std::vector<int> sa;
void work(){
    cin >> str;
    str += '\0';
    SAIS::sais(str, str.size(), 128, sa);
    int len = sa.size();
    UP(i, 1, len){
        cout << sa[i]+1 << ' ';
    }
    cout << '\n';
}
} // {}}}
int main(){m::work(); return 0;}

性能测试

使用的文本为 command_block 大佬的《多项式计数杂谈》(约 128KiB),复读若干次。字符集大小为 65536。

对照组为某谷 P3809 最高赞题解,双方均用 fread 读入后转 wchar_t,用 printf 输出,开 -O3

用时的单位为秒。

未经多次测试,结果可能不准确。

文件大小 SA-IS 倍增
100K 0.016 0.031
1.5M 0.173 0.359
20M 3.369 21.271

推荐阅读

A walk through the SA-IS Suffix Array Construction Algorithm(好文,从这学的)