usestd::cmp::max;usestd::cmp::Ordering;usestd::slice::from_raw_parts_mut;constLTYPE:bool=false;constSTYPE:bool=true;constMAX_SA_VALUE:usize=usize::MAX/2;constEMPTY:usize=MAX_SA_VALUE+1;constUNIQUE:usize=MAX_SA_VALUE+2;constMULTI:usize=MAX_SA_VALUE+3;// >= 258fnlms_str_cmp<E:Ord>(l1:&[E],l2:&[E])->Ordering{for(x,y)inl1.iter().zip(l2.iter()){letcmp_res=x.cmp(&y);ifcmp_res!=Ordering::Equal{returncmp_res;}}Ordering::Equal}#[inline]fnpat_char_type(cur:usize,prev:usize,last_scanned_type:bool)->bool{ifcur<prev||cur==prev&&last_scanned_type==STYPE{STYPE}else{LTYPE}}fnrename_pat(pat:&mut[usize],sa:&mut[usize]){letpatlastpos=pat.len()-1;// 全部刷成bucket head//sa.fill(0);foriin0..sa.len(){sa[i]=0}foriin0..pat.len(){sa[pat[i]]+=1}foriin1..sa.len(){sa[i]+=sa[i-1]}foriin0..pat.len()-1{pat[i]=sa[pat[i]]-1;};// 将L-suffix刷成bucket head//sa.fill(0);foriin0..sa.len(){sa[i]=0}foriin0..pat.len(){sa[pat[i]]+=1}letmutlast_scanned_type=STYPE;pat[patlastpos]=0;foriin(0..pat.len()-1).rev(){ifpat_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;}}}fnsort_lms_char(pat:&mut[usize],sa:&mut[usize])->usize{//sa.fill(EMPTY);foriin0..sa.len(){sa[i]=EMPTY}letmutlast_scanned_type=STYPE;foriin(0..pat.len()-1).rev(){ifpat_char_type(pat[i],pat[i+1],last_scanned_type)==STYPE{last_scanned_type=STYPE;}else{iflast_scanned_type==STYPE{// pat[i + 1] is LMS typesa[pat[i+1]]+=1;}last_scanned_type=LTYPE;}}letmutlms_cnt=0;last_scanned_type=STYPE;foriin(0..pat.len()-1).rev(){ifpat_char_type(pat[i],pat[i+1],last_scanned_type)==STYPE{last_scanned_type=STYPE;}else{lete_i=i+1;lete=pat[e_i];iflast_scanned_type==STYPE{// pat[i + 1] is LMS typelms_cnt+=1;ifsa[e]==UNIQUE{sa[e]=e_i;}elseifsa[e]>=MULTI&&sa[e-1]==EMPTY{ifsa[e-2]==EMPTY{sa[e-2]=e_i;sa[e-1]=1;// set counter}else{// MUL = 2sa[e]=e_i;sa[e-1]=EMPTY;}}elseifsa[e]>=MULTI&&sa[e-1]!=EMPTY{letc=sa[e-1];// get counterifsa[e-2-c]==EMPTY{sa[e-2-c]=e_i;sa[e-1]+=1;// update counter}else{forjin(1..c+1).rev(){sa[e-c+j]=sa[e-2-c+j]}sa[e-c]=e_i;sa[e-c-1]=EMPTY;}}elseifsa[e]<EMPTY{forjin(0..e).rev(){ifsa[j]==EMPTY{sa[j]=e_i;break;}}}}last_scanned_type=LTYPE;}}foriin(0..pat.len()).rev(){ifsa[i]>=MULTI{letc=sa[i-1];forjin(1..c+1).rev(){// 逆序防止前面的覆盖后面的sa[i-c+j]=sa[i-2-c+j];}sa[i-c-1]=EMPTY;sa[i-c]=EMPTY;}}lms_cnt}fnsort_lms_substr(pat:&mut[usize],sa:&mut[usize]){// step 1induced_sort(pat,sa);// step 2letpat_last_pos=pat.len()-1;letmutlms_cnt=0;letmuti=pat_last_pos;letmutbucket_tail_ptr=pat_last_pos+1;// for renamed bucket verletmutbucket=EMPTY;// 可以省略,但是为了书写代码方便letmutnum=0;// S type number of bucketwhilei>0{ifpat[sa[i]]!=bucket{// reach new bucketnum=0;letmutl=0;whilepat[sa[i-l]]==pat[sa[i]]{// 扫描桶来计算桶中S字符数量,根据定义 当l=i时循环必然终止letpat_i=sa[i-l];// l < i, 即 i - l > 0, 0 <= pat_i < patlen - 1ifpat[pat_i]<pat[pat_i+1]{letmutk=pat_i;whilek>0&&pat[k-1]==pat[pat_i]{k-=1}num+=pat_i-k+1;}else{break;// bucket不含S字符,结束扫描}l+=1;}bucket_tail_ptr=i;bucket=pat[sa[bucket_tail_ptr]];}ifnum>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 = 0lms_cnt+=1;//sa[0..pat_last_pos - lms_cnt + 1].fill(EMPTY);foriin0..pat_last_pos-lms_cnt+1{sa[i]=EMPTY}}fnconstruct_pat1(pat:&mut[usize],sa:&mut[usize],lms_cnt:usize)->bool{letpatlen=pat.len();letmutprev_lms_str_len=1;letmutrank=0;sa[(patlen-1)/2]=rank;letmuthas_duplicated_char=false;foriinpatlen-lms_cnt+1..patlen{// 从警戒哨字符的下一个字符开始letmutj=sa[i];whilepat[j]<=pat[j+1]{j+=1}// 寻找suf(sa[i])右边第一个L字符,因为排除了警戒哨这个LMS后缀,所以必然不会越界letmutk=j;whilek+1<patlen&&pat[k]>=pat[k+1]{k+=1}// 找到suf(sa[i])右边第一个LMS字符letcur_lms_str_len=k+1-sa[i];letcmp_res=lms_str_cmp(&pat[sa[i]..sa[i]+cur_lms_str_len],&pat[sa[i-1]..sa[i-1]+prev_lms_str_len]);ifcmp_res!=Ordering::Equal{rank+=1}ifrank==sa[sa[i-1]/2]{has_duplicated_char=true;}letrank_index=sa[i]/2;sa[rank_index]=rank;// 整除prev_lms_str_len=cur_lms_str_len;}// move to head of saletmutj=0;foriin0..patlen-lms_cnt{ifsa[i]!=EMPTY{sa[j]=sa[i];ifi>j{sa[i]=EMPTY;}j+=1;}}//sa[lms_cnt..patlen].fill(EMPTY);foriinlms_cnt..patlen{sa[i]=EMPTY}has_duplicated_char}fnsort_lms_suf(pat:&mut[usize],sa:&mut[usize],lms_cnt:usize,has_duplicated_char:bool){// solve T1 recursivelyletpatlen=pat.len();letsalen=sa.len();unsafe{letsa_ptr=sa.as_mut_ptr();letmutpat1=from_raw_parts_mut(sa_ptr,lms_cnt);letmutsa1=from_raw_parts_mut(sa_ptr.offset((patlen-lms_cnt)asisize),salen-(patlen-lms_cnt));ifhas_duplicated_char{_compute_suffix_array_16_1(&mutpat1,&mutsa1);}else{foriin0..lms_cnt{sa1[pat1[i]]=i}}}// move SA1 to SA[0...n1-1]foriin0..lms_cnt{sa[i]=sa[patlen-lms_cnt+i];}// put all LMS-suffixes in SA tailletmutlast_scanned_type=STYPE;letmutj=0;foriin(0..pat.len()-1).rev(){ifpat[i]<pat[i+1]||pat[i]==pat[i+1]&&last_scanned_type==STYPE{last_scanned_type=STYPE;}else{iflast_scanned_type==STYPE{sa[patlen-1-j]=i+1;j+=1;}last_scanned_type=LTYPE;}}// backward map the LMS-suffixes rankforiin0..lms_cnt{letrelative_rank=sa[i];sa[i]=sa[patlen-lms_cnt+relative_rank];sa[patlen-lms_cnt+relative_rank]=EMPTY;}letmuttail=EMPTY;letmutrfp=EMPTY;foriin(1..lms_cnt).rev(){// sa[0] 保持原位ifpat[sa[i]]!=tail{tail=pat[sa[i]];rfp=tail;}sa[rfp]=sa[i];ifrfp!=i{sa[i]=EMPTY}rfp-=1;}}// PASS!fninduced_sort(pat:&mut[usize],sa:&mut[usize]){letpatlen=pat.len();// place L-suff in SA// initletmutlast_scanned_type=STYPE;foriin(0..patlen-1).rev(){ifpat_char_type(pat[i],pat[i+1],last_scanned_type)==LTYPE{sa[pat[i]]+=1;// >= EMPTYlast_scanned_type=LTYPE;}else{last_scanned_type=STYPE;}}//placeletmuti=0;whilei<patlen{ifsa[i]<EMPTY&&sa[i]>0{letj=sa[i]-1;letmutis_ltype=false;ifpat[j]>pat[j+1]{is_ltype=true;}elseifpat[j]==pat[j+1]{// 判断sa[i]是否是L后缀的编号letnext_i=sa[pat[sa[i]]];ifnext_i>=MULTI{is_ltype=true;}elseifnext_i<EMPTY&&pat[sa[i]]+1<patlen{ifsa[pat[sa[i]]+1]==EMPTY{is_ltype=true;}elseifsa[pat[sa[i]]+1]<EMPTY{ifpat[sa[pat[sa[i]]+1]]==pat[sa[i]]{is_ltype=true;}}}}ifis_ltype{ifsa[pat[j]]==UNIQUE{sa[pat[j]]=j;}elseifsa[pat[j]]>=MULTI&&sa[pat[j]+1]==EMPTY{ifsa[pat[j]]-EMPTY>2{sa[pat[j]+2]=j;sa[pat[j]+1]=1;// set counter}else{sa[pat[j]]=j;}}elseifsa[pat[j]]>=MULTI&&sa[pat[j]+1]!=EMPTY{lete=pat[j];letc=sa[e+1];letlfp=e+c+2;ifc+2<sa[pat[j]]-EMPTY{// 没到bucket尾部sa[lfp]=j;sa[e+1]+=1;// update counter}else{forkin1..c+1{sa[e+k-1]=sa[e+k+1];}sa[e+c]=j;sa[e+c+1]=EMPTY;ifi>=e+2&&i<=e+c+1{i-=2;}}}elseifsa[pat[j]]<EMPTY{forkinpat[j]..patlen{ifsa[k]==EMPTY{sa[k]=j;break;}}}}}elseifsa[i]>=MULTI{i+=1;}i+=1;}// remove LMS-suff form SA, 一个桶里可能有多个LMS后缀last_scanned_type=STYPE;foriin(0..pat.len()-1).rev(){ifpat_char_type(pat[i],pat[i+1],last_scanned_type)==STYPE{last_scanned_type=STYPE;}else{iflast_scanned_type==STYPE{// pat[i + 1] is LMS typeifsa[pat[i+1]]<=EMPTY{sa[pat[i+1]]=UNIQUE;}else{sa[pat[i+1]]+=1;}}last_scanned_type=LTYPE;}}i=patlen-1;whilei>0{ifsa[i]>EMPTY{letc=sa[i]-EMPTY;forkin0..c{sa[i-k]=EMPTY;}i-=c-1;}i-=1;}sa[0]=pat.len()-1;// place S-suff in SA// initletmutlast_scanned_type=STYPE;foriin(0..patlen-1).rev(){ifpat_char_type(pat[i],pat[i+1],last_scanned_type)==STYPE{ifsa[pat[i]]>=EMPTY{sa[pat[i]]+=1;}else{sa[pat[i]]=UNIQUE;}last_scanned_type=STYPE;}else{last_scanned_type=LTYPE;}}i=patlen-1;whilei>0{ifsa[i]<EMPTY&&sa[i]>0{letj=sa[i]-1;letmutis_stype=false;ifpat[j]<pat[j+1]{is_stype=true;}elseifpat[j]==pat[j+1]{// 判断sa[i]是否是S后缀的编号letnext_i=sa[pat[sa[i]]];ifnext_i>=MULTI{is_stype=true;}elseifnext_i<EMPTY&&pat[sa[i]]-1>0{ifsa[pat[sa[i]]-1]==EMPTY{is_stype=true;}elseifsa[pat[sa[i]]-1]<EMPTY{ifpat[sa[pat[sa[i]]-1]]==pat[sa[i]]{is_stype=true;}}}}ifis_stype{ifsa[pat[j]]==UNIQUE{sa[pat[j]]=j;}elseifsa[pat[j]]>=MULTI&&sa[pat[j]-1]==EMPTY{ifsa[pat[j]]-EMPTY>2{sa[pat[j]-2]=j;sa[pat[j]-1]=1;// set counter}else{sa[pat[j]]=j;}}elseifsa[pat[j]]>=MULTI&&sa[pat[j]-1]!=EMPTY{lete=pat[j];letc=sa[e-1];letnum=sa[pat[j]]-EMPTY;ifc+2<num{// 没到bucket头部letrfp=e-c-2;sa[rfp]=j;sa[e-1]+=1;}else{forkin1..c+1{sa[e-k+1]=sa[e-k-1];}sa[e-c]=j;sa[e-c-1]=EMPTY;ifi>=e-num+1&&i<=e-2{i+=2;}}}elseifsa[pat[j]]<EMPTY{forkin(0..pat[j]).rev(){ifsa[k]==EMPTY{sa[k]=j;break;}}}}}elseifsa[i]>=MULTI{i-=1;}i-=1;}}fn_compute_suffix_array_16_1(pat:&mut[usize],sa:&mut[usize]){rename_pat(pat,sa);letlms_cnt=sort_lms_char(pat,sa);sort_lms_substr(pat,sa);lethas_duplicated_char=construct_pat1(pat,sa,lms_cnt);sort_lms_suf(pat,sa,lms_cnt,has_duplicated_char);induced_sort(pat,sa);}pubfnsuffix_array_16(pat:&[u8])->Vec<usize>{letmutpat=pat.into_iter().map(|x|*xasusize).collect::<Vec<usize>>();pat.push(0);letmutsa=vec![0;max(pat.len(),256)*1];_compute_suffix_array_16_1(&mutpat[..],&mutsa[..]);sa}fninput()->String{usestd::io;letmutinput=String::new();io::stdin().read_line(&mutinput).unwrap();String::from(input.trim())}fnmain(){letpat=input();letsa_16=suffix_array_16(pat.as_bytes());foriin1..pat.len()+1{print!("{} ",sa_16[i]+1)}}
在只读的整形字母表上的后缀排序
使用复杂方法解决复杂问题,通过分治,解决空间紧张的问题。
算法实现的难点在于在 上构建 BitMaps5,来替代本来由重命名后的 T 所指示的指示桶尾/桶头的位置。
这里的 BitMaps 指得是使用比特向量(bit vector)表示的有序字典(multiset),是一种紧凑型结构(compact data structure)。
有兴趣了解的暂时只能阅读原文以及本文引用的 BitMaps 的有关论文自行了解。
在只读的一般字母表上的后缀排序
前置知识是归并排序和堆排序。
由于笔者对于其中确定字符类型的方法的时间复杂度有疑问,这里也不再介绍,建议阅读原文自行了解。
注解
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. 11147. Springer. pp. 268–284. arXiv:1610.08305. doi:10.1007/978-3-030-00479-8_22. ISBN:978-3-030-00478-1. ↩↩
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. ↩
如果是 LML 后缀,就先诱导 S 型后缀,唯一区别是计算 LML 后缀时需要将警戒哨也算进去。 ↩
Gonzalo Navarro and Eliana Providel. Fast, small, simple rank/select on bitmaps. In Proc. 11th International Symposium on Experimental Algorithms (SEA), pages 295–306, 2012. ↩