文字列検索アルゴリズム入門 - KMP・Boyer-Moore・Rabin-Karpを動くコードで理解する

文字列検索アルゴリズム入門 - KMP・Boyer-Moore・Rabin-Karpを動くコードで理解する

作成日:
読了:41
更新日:

エディタの検索ボックスに文字を打ち込むと、数MBのファイルからでも一瞬で該当箇所がハイライトされます。grep は数GBのログを秒単位でなめ切ります。この「文字列の中から文字列を探す」という、あまりに日常的すぎて意識すらしない処理の裏側には、50年近く磨かれてきたアルゴリズムの蓄積があります。この記事では、素朴な全探索から始めて、Rabin-Karp(ローリングハッシュ)、KMPBoyer-MooreZ-algorithm という4つの古典を、実際に動くPythonコードで順に組み立てていきます。最後に、CPythonの str.find やgrep / ripgrepが実際に何をしているのかという実務の話まで踏み込みます。二分探索入門Union-Find入門と同じアルゴリズム入門シリーズの1本です。

問題設定と用語 — テキストからパターンを探す

扱う問題は厳密文字列照合(exact string matching)と呼ばれるものです。

  • テキスト T: 探される側の長い文字列。長さを n とします
  • パターン P: 探す側の短い文字列。長さを m とします
  • 求めたいもの: T の中で P と完全に一致する出現位置のすべて

「位置」は T の何文字目から P が始まるかを0始まりのインデックスで表します。たとえば TabracadabraPabra なら、答えは [0, 7] です。

いくつか約束を決めておきます。

  • 出現は重なってよい: TaaaaaaPaa のとき、答えは [0, 1, 2, 3, 4] の5個です。重なりを許さない数え方(Pythonの str.count はこちら)とは結果が変わります
  • 最初の1個だけでよい場合もある: 実務では「あるか無いか」「最初の位置はどこか」で十分なことが多く、この場合は見つかった時点で打ち切れます
  • 比較単位はここでは文字: 実際のライブラリはバイト単位で動くものが多く、マルチバイト文字の扱いが別問題として出てきます(詳しくは文字コードとUTF-8入門を参照してください)

なお本記事では、正規表現のようなパターン言語ではなく、あくまで固定文字列の照合を扱います。正規表現との住み分けは記事の後半で整理します。

素朴な全探索 — まずはここから

もっとも素直な方法は、T の各位置に P を置いてみて、先頭から1文字ずつ比べることです。

def brute_force_search(text, pattern):
    n, m = len(text), len(pattern)
    result = []
    for s in range(n - m + 1):          # 開始位置 s を左から順に試す
        j = 0
        while j < m and text[s + j] == pattern[j]:   # 前から1文字ずつ照合
            j += 1
        if j == m:                      # 最後まで一致したら出現位置として記録
            result.append(s)
    return result
 
 
print(brute_force_search("abracadabra", "abra"))
print(brute_force_search("aaaaaa", "aa"))
print(brute_force_search("hello", "xyz"))
[0, 7]
[0, 1, 2, 3, 4]
[]

正しく動きますが、計算量は最悪 O(n * m) です。開始位置が n - m + 1 通りあり、それぞれで最大 m 回比較するからです。

問題は「最悪ケースが本当に起きるのか」です。英文のような普通のテキストでは、たいてい1〜2文字目で不一致になって次へ進むので、実測はほぼ O(n) に見えます。しかし、同じ文字が延々と続くテキストでは話が変わります。

def count_comparisons(text, pattern):
    n, m = len(text), len(pattern)
    comparisons = 0
    for s in range(n - m + 1):
        j = 0
        while j < m:
            comparisons += 1                       # 1文字比較するたびに数える
            if text[s + j] != pattern[j]:
                break
            j += 1
    return comparisons
 
 
print(count_comparisons("a" * 20, "aaab"))         # 最悪ケース
print(count_comparisons("abcdefghijklmnopqrst", "aaab"))   # 普通のケース
print(count_comparisons("a" * 1000, "a" * 9 + "b"))
68
18
9910

aaaa...a から aaab を探すと、どの開始位置でも「3文字一致して4文字目で不一致」となり、毎回 m 回近い比較を払わされます。20文字のテキストで68回、1000文字で9910回とほぼ n * m に張り付いています。一方、普通のテキストでは18回、つまりほぼ n 回です。

この「同じ文字の繰り返し」はDNA配列、バイナリデータ、区切り文字が連続するログなど、現実に十分あり得ます。この最悪ケースを潰すのが以降のアルゴリズムの目的です。計算量の記法そのものが不安な方は、先に計算量とBig-O記法 入門に目を通しておくと読みやすくなります。

Rabin-Karp — ハッシュで窓を滑らせる

Rabin-Karpアルゴリズムは、Richard M. KarpとMichael O. Rabinが1987年に発表した手法です("Efficient randomized pattern-matching algorithms", IBM Journal of Research and Development 31(2), pp.249-260)。発想は「文字列どうしを比べる前に、ハッシュ値どうしを比べてほとんどの位置を門前払いする」というものです。

ローリングハッシュの原理

素直にやると、各位置で長さ m の部分文字列のハッシュを計算し直すことになり、結局 O(n * m) です。ここで効くのがローリングハッシュという考え方です。

文字列を基数 B の多項式と見なします。長さ m の窓 T[s..s+m-1] のハッシュを次のように定義します。

H(s) = T[s]*B^(m-1) + T[s+1]*B^(m-2) + ... + T[s+m-1]*B^0   (mod M)

すると、窓を1つ右にずらしたときのハッシュは、先頭の文字を引き、全体をB倍し、新しい末尾の文字を足すだけで求まります。

H(s+1) = (H(s) - T[s]*B^(m-1)) * B + T[s+m]   (mod M)

つまり窓の更新が O(1) です。ハッシュそのものの考え方はハッシュテーブル入門と共通ですが、こちらは「隣の入力へ安く移れる」性質を積極的に使う点が特徴です。

実装

ハッシュが一致しても本当に一致しているとは限らないので、必ず実文字列で検証します。これを省くと誤検出が起きます。

def rabin_karp_search(text, pattern, base=256, mod=(1 << 61) - 1):
    n, m = len(text), len(pattern)
    if m == 0 or n < m:
        return []
    high = pow(base, m - 1, mod)        # 先頭文字にかかる係数 B^(m-1)
    ph = th = 0
    for i in range(m):                  # パターンと最初の窓のハッシュを計算
        ph = (ph * base + ord(pattern[i])) % mod
        th = (th * base + ord(text[i])) % mod
    result = []
    for s in range(n - m + 1):
        if th == ph and text[s:s + m] == pattern:   # ハッシュ一致は候補どまり。実文字列で検証する
            result.append(s)
        if s + m < n:                   # 窓を1つ右へずらす(O(1) で更新)
            th = (th - ord(text[s]) * high) % mod
            th = (th * base + ord(text[s + m])) % mod
    return result
 
 
print(rabin_karp_search("abracadabra", "abra"))
print(rabin_karp_search("aaaaaa", "aa"))
print(rabin_karp_search("the quick brown fox", "quick"))
[0, 7]
[0, 1, 2, 3, 4]
[4]

計算量と、基数・modの選び方

前処理は O(m)、検索は平均 O(n)最悪 O(n * m) です。追加メモリは O(1) で済みます。最悪ケースはハッシュ衝突が毎回起きる場合で、そのとき検証コストが全位置で発生します。

したがって実用性はハッシュ関数の質にすべて依存します。押さえるべき点は次の通りです。

  • mod は大きな素数を選ぶ: 上の実装では 2^61 - 1(メルセンヌ素数)を使っています。32ビット程度の小さな法だと、大きな入力では偶然の衝突が現実的な頻度で起きます
  • 基数は文字集合の大きさ以上に取る: バイト列なら256以上。理論的な保証を重視するなら、実行のたびに基数を乱数で選ぶ(ランダム化する)と、特定の入力を狙い撃ちする攻撃を防げます
  • ハッシュ一致を「一致」と即断しない: 検証を省いた実装は、悪意ある入力に対して誤った結果を返します。ハッシュ値だけで判定する設計は避けてください

複数パターンへの拡張

Rabin-Karpの真価は長さの揃った複数パターンを一度に探すときに出ます。パターン側のハッシュを集合に入れておき、窓のハッシュがその集合に含まれるかを見るだけです。

def rabin_karp_multi(text, patterns, base=256, mod=(1 << 61) - 1):
    m = len(patterns[0])
    assert all(len(p) == m for p in patterns)      # 長さが揃っている前提
    table = {}
    for p in patterns:                             # ハッシュ値 -> 同じハッシュを持つパターン一覧
        h = 0
        for ch in p:
            h = (h * base + ord(ch)) % mod
        table.setdefault(h, []).append(p)
    n = len(text)
    if n < m:
        return []
    high = pow(base, m - 1, mod)
    th = 0
    for i in range(m):
        th = (th * base + ord(text[i])) % mod
    result = []
    for s in range(n - m + 1):
        if th in table:                            # 候補が出たときだけ実文字列で確認
            window = text[s:s + m]
            for p in table[th]:
                if window == p:
                    result.append((s, p))
        if s + m < n:
            th = (th - ord(text[s]) * high) % mod
            th = (th * base + ord(text[s + m])) % mod
    return result
 
 
print(rabin_karp_multi("the cat sat on the mat", ["cat", "mat", "hat"]))
[(4, 'cat'), (19, 'mat')]

k個のパターンを個別に探すと O((n + m) * k) かかるところ、この方式なら期待 O(n + k * m) です。剽窃検出のように「大量の短い断片を一度に照合する」用途で使われる理由がここにあります。

KMP — ずらし幅をパターンから前計算する

Knuth-Morris-Prattアルゴリズム(KMP)は、Donald Knuth、James H. Morris、Vaughan Prattの3名が1977年に発表しました("Fast pattern matching in strings", SIAM Journal on Computing 6(2), pp.323-350)。Wikipediaによれば、Morrisが先に考案し、Knuthがオートマトン理論から数週間後に独立に再発見したという経緯です。

発想: 不一致から情報を引き出す

素朴な全探索の無駄はどこにあったか。aaaa に対して aaab を照合し、4文字目で失敗したとします。このとき私たちは「直前の3文字が aaa だった」ことをすでに知っているのに、開始位置を1つ進めてまた先頭から比べ直しています。

KMPは、この「すでに一致していた接頭辞」の情報を捨てません。そのために、パターンだけを見て前計算する配列を用意します。

LPS配列(failure function)

lps[i] を、P[0..i] の接頭辞であり同時に接尾辞でもある文字列のうち、最長のもの(ただし全体そのものを除く)の長さと定義します。LPSはLongest Proper Prefix which is also Suffixの略で、教科書によってはfailure functionやprefix functionとも呼ばれます。

def build_lps(pattern):
    m = len(pattern)
    lps = [0] * m
    k = 0                                  # 直前の位置までで一致している接頭辞の長さ
    for i in range(1, m):
        while k > 0 and pattern[i] != pattern[k]:
            k = lps[k - 1]                 # 合わなければ1段短い候補へ後退する
        if pattern[i] == pattern[k]:
            k += 1
        lps[i] = k
    return lps
 
 
print(build_lps("ababaca"))
print(build_lps("aabaaab"))
print(build_lps("abcdef"))
[0, 0, 1, 2, 3, 0, 1]
[0, 1, 0, 1, 2, 2, 3]
[0, 0, 0, 0, 0, 0]

ababaca の5文字目まで(ababa)を見ると、接頭辞 aba と接尾辞 aba が一致するので lps[4] = 3 です。すべての文字が異なる abcdef では、どこにも一致がないので全部0になります。

検索本体

LPS配列があれば、不一致が起きたときにパターン側のポインタだけを後退させ、テキスト側のポインタは一度も戻しません

def kmp_search(text, pattern):
    n, m = len(text), len(pattern)
    if m == 0:
        return list(range(n + 1))
    lps = build_lps(pattern)
    result = []
    k = 0                                  # 現在一致しているパターンの長さ
    for i in range(n):                     # テキストのポインタ i は決して戻らない
        while k > 0 and text[i] != pattern[k]:
            k = lps[k - 1]                 # 不一致ならパターン側だけ後退
        if text[i] == pattern[k]:
            k += 1
        if k == m:                         # 全部一致した
            result.append(i - m + 1)
            k = lps[k - 1]                 # 重なる出現も拾えるように後退させる
    return result
 
 
print(kmp_search("abracadabra", "abra"))
print(kmp_search("aaaaaa", "aa"))
print(kmp_search("ababababa", "ababa"))
[0, 7]
[0, 1, 2, 3, 4]
[0, 2, 4]

ababababa から ababa を探すと、位置0、2、4で重なりながら3回見つかります。一致直後に k = lps[k-1] としているので、重なる出現も取りこぼしません。

なぜ O(n + m) になるのか

計算量は前処理 Θ(m)、検索 Θ(n)、追加メモリ O(m) です。素朴な全探索と違い、最悪ケースでもこの上界を超えません。

理由はならし計算量の議論で説明できます。検索ループでは i が n 回進み、その各回で k はたかだか1しか増えません。一方 while 内の k = lps[k-1]k を必ず減らします。k は増える総量が高々 n なので、減る総量も高々 n です。したがって while の総実行回数は n 回を超えられず、全体で O(n) に収まります。LPS配列の構築も同じ議論で O(m) です。

「テキストのポインタが戻らない」という性質は、ファイルやネットワークからのストリームをそのまま流し込めることを意味します。バッファに全体を溜めなくても照合できるのは、KMPの実装上の大きな利点です。

Boyer-Moore — 後ろから比較して大きく飛ぶ

Boyer-Mooreアルゴリズムは、Robert S. BoyerとJ Strother Mooreが1977年に発表しました("A fast string searching algorithm", Communications of the ACM 20(10), pp.762-772)。KMPと同じ年ですが、発想はほぼ正反対です。

Boyer-Mooreはパターンの後ろから比較します。この一点で何が変わるかというと、テキストの大半の文字を一度も見ずに飛ばせるようになります。

bad character rule

窓の中でパターンの位置 j が不一致になり、そこにあったテキストの文字が c だったとします。このとき、パターンの中で c が最後に現れる位置を調べます。

  • パターンに c が含まれるなら、その位置が j に重なるように窓をずらす
  • パターンに cまったく含まれないなら、窓全体を c の右へ飛ばしてよい

後者が強力です。パターンに無い文字を1つ踏むだけで、m 文字まるごとスキップできます。

good suffix rule

もう1つの規則がgood suffix ruleです。不一致が起きた時点で、すでに末尾の何文字かは一致しています。その一致した接尾辞が、パターンの別の場所にも現れているなら、そこが重なるようにずらします。どこにも現れないなら、「一致した接尾辞の接尾辞」がパターンの接頭辞になっている最長のところまで飛ばせます。

実装は2つの前処理表からなります。前処理は Θ(m) 時間、空間は Θ(k + m)(k はアルファベットサイズ)です。

def build_good_suffix(pattern):
    m = len(pattern)
    shift = [0] * (m + 1)
    border = [0] * (m + 1)
    i, j = m, m + 1
    border[i] = j
    while i > 0:                            # ケース1: 一致した接尾辞が他の場所にも現れる
        while j <= m and pattern[i - 1] != pattern[j - 1]:
            if shift[j] == 0:
                shift[j] = j - i
            j = border[j]
        i -= 1
        j -= 1
        border[i] = j
    j = border[0]
    for i in range(m + 1):                  # ケース2: 接尾辞の一部がパターンの接頭辞になる
        if shift[i] == 0:
            shift[i] = j
        if i == j:
            j = border[j]
    return shift
 
 
def boyer_moore_search(text, pattern):
    n, m = len(text), len(pattern)
    if m == 0 or n < m:
        return []
    last = {}
    for i, ch in enumerate(pattern):
        last[ch] = i                        # bad character rule 用: 各文字が最後に現れる位置
    gs = build_good_suffix(pattern)
    result = []
    s = 0
    while s <= n - m:
        j = m - 1
        while j >= 0 and pattern[j] == text[s + j]:   # 後ろから前へ比較する
            j -= 1
        if j < 0:
            result.append(s)
            s += gs[0]
        else:
            bad_char = j - last.get(text[s + j], -1)  # 未登場の文字なら -1 として大きく飛ぶ
            s += max(gs[j + 1], bad_char)             # 2つの規則のうち大きいほうを採用
    return result
 
 
print(build_good_suffix("ANPANMAN"))
print(boyer_moore_search("abracadabra", "abra"))
print(boyer_moore_search("HERE IS A SIMPLE EXAMPLE", "EXAMPLE"))
[6, 6, 6, 6, 6, 6, 3, 8, 1]
[0, 7]
[17]

Horspool版 — good suffix ruleを捨てる

good suffix ruleは正しく書くのが難しく、実装バグの温床でもあります。そこでNigel Horspoolが1980年に提案したのが、bad character ruleだけを使う簡略版です("Practical fast searching in strings", Software: Practice and Experience 10(6), pp.501-506)。

Horspool版の工夫は、ずらし幅を不一致が起きた位置の文字ではなく、窓の末尾の文字で引くことです。これによりテーブルが1本で済み、コードが劇的に短くなります。

def horspool_search(text, pattern):
    n, m = len(text), len(pattern)
    if m == 0 or n < m:
        return []
    shift = {}
    for i in range(m - 1):                  # 最後の1文字は入れない(入れると0になってしまう)
        shift[pattern[i]] = m - 1 - i
    result = []
    s = 0
    while s <= n - m:
        j = m - 1
        while j >= 0 and pattern[j] == text[s + j]:
            j -= 1
        if j < 0:
            result.append(s)
        s += shift.get(text[s + m - 1], m)  # 窓の末尾の文字でずらし幅を引く。未登場なら m 飛ばす
    return result
 
 
print(horspool_search("abracadabra", "abra"))
print(horspool_search("aaaaaa", "aa"))
print(horspool_search("HERE IS A SIMPLE EXAMPLE", "EXAMPLE"))
[0, 7]
[0, 1, 2, 3, 4]
[17]

EXAMPLE のずらし表は A:4, E:6, L:1, M:3, P:2, X:5 となります。この表で24文字のテキストを検索すると、調べた窓は5個、文字比較は15回で済みます。素朴な全探索は27回の比較(窓は18個)でした。テキストが長くなるほどこの差は開きます。

なぜ実務で速く、どこで壊れるか

Boyer-Mooreの最良ケースは Ω(n / m) 比較です。テキストの文字数より少ない回数で終わる、つまり劣線形(sublinear)に振る舞い得るという点が、他のアルゴリズムにない特徴です。パターンが長いほど、そしてアルファベットが大きい(日本語や自然文のように文字種が多い)ほど有利になります。

一方で弱点もはっきりしています。

  • Galil ruleを入れないと最悪 O(n * m): パターンがテキスト中に何度も出現する場合に劣化します。Galil ruleを組み込むと全ケースで線形 O(n + m) になります
  • アルファベットが小さいと効かない: DNA配列(4文字)やビット列では、bad character ruleがほとんどスキップを生みません
  • Horspool版の最悪は O(n * m): 平均・最良は Θ(n) ですが、最悪の保証はありません

Z-algorithm — Z配列ひとつで検索する

Z-algorithmは、KMPとは別の角度から線形時間照合を実現する手法です。Dan Gusfieldの教科書『Algorithms on Strings, Trees, and Sequences』(1997年)で "fundamental preprocessing" として解説され、競技プログラミングの世界で広く使われるようになりました。なお、この手法の初出をMain and Lorentz(1984年, Journal of Algorithms 5(3), pp.422-432)に帰す記述をよく見かけますが、筆者が確認できた範囲では一次情報にたどり着けなかったため未確認としておきます。

Z配列の定義

文字列 S に対して、z[i]SS[i..] の最長共通接頭辞の長さと定義します。z[0] は慣習として len(S) とします。

線形構築の肝は、直近に見つかった「Sの接頭辞と一致する区間」を使い回すことです。この区間を [left, right) として保持し、新しい位置 i がその中に入っていれば、対称な位置 i - left のZ値をそのまま初期値として流用できます。

def z_array(s):
    n = len(s)
    z = [0] * n
    if n == 0:
        return z
    z[0] = n
    left = right = 0                        # 現在わかっている「接頭辞と一致する区間」
    for i in range(1, n):
        if i < right:
            z[i] = min(right - i, z[i - left])   # 対称な位置の結果を流用して比較を省く
        while i + z[i] < n and s[z[i]] == s[i + z[i]]:
            z[i] += 1                       # はみ出した分だけ素朴に伸ばす
        if i + z[i] > right:
            left, right = i, i + z[i]       # 区間を更新
    return z
 
 
print(z_array("aabxaayaab"))
print(z_array("aaaaa"))
print(z_array("aabcaabxaaaz"))
[10, 1, 0, 0, 2, 1, 0, 3, 1, 0]
[5, 4, 3, 2, 1]
[12, 1, 0, 0, 3, 1, 0, 0, 2, 2, 1, 0]

while ループで right は単調増加し、全体を通して n を超えないので、構築は O(n) 時間・O(n) 空間です。

パターン検索への応用

Z配列そのものはパターン検索の道具に見えませんが、パターンとテキストを区切り文字でつないで1本の文字列にするという一手で検索になります。

def z_search(text, pattern, sep="\x00"):
    assert sep not in text and sep not in pattern   # 区切り文字は両方に出現しないもの
    m = len(pattern)
    z = z_array(pattern + sep + text)               # 「パターン + 区切り + テキスト」を作る
    return [i - m - 1 for i in range(m + 1, len(z)) if z[i] == m]  # z値がmならそこが出現位置
 
 
print(z_search("abracadabra", "abra"))
print(z_search("ababababa", "ababa"))
[0, 7]
[0, 2, 4]

区切り文字が両方に現れないことが必須条件です。ここを守らないと、z[i] が m を超えてしまい判定が壊れます。テキストが任意のバイト列を含み得る場面では、区切り文字を安全に選べないことがある点に注意してください。

計算量は前処理と検索をまとめて O(n + m)、追加メモリは O(n + m)。KMPと同じ線形ですが、連結した文字列ぶんのメモリを使うのが違いです。その代わり、Z配列は「各位置から見た接頭辞との一致長」という汎用的な情報なので、周期性の検出や文字列の圧縮表現など、検索以外の応用にもそのまま使えます。

5つのアルゴリズムの比較

ここまでの内容を1枚にまとめます。n はテキスト長、m はパターン長、k はアルファベットサイズです。

アルゴリズム前処理検索(平均)検索(最悪)追加メモリ得意な場面
素朴な全探索なしほぼ O(n)O(n * m)O(1)短いテキスト、実装を最優先したいとき
Rabin-KarpO(m)O(n)O(n * m)O(1)同じ長さの複数パターン、剽窃検出
KMPΘ(m)Θ(n)Θ(n)O(m)最悪保証が要る、ストリーム処理
Boyer-MooreΘ(m)劣線形になり得るO(n * m)(Galil ruleで O(n + m)Θ(k + m)長いパターン、文字種が多いテキスト
Boyer-Moore-HorspoolO(m)Θ(n)O(n * m)Θ(k)実装を簡単にしたいとき、一般的なテキスト
Z-algorithm連結して O(n + m)O(n + m)O(n + m)O(n + m)前処理情報を他の用途にも使いたいとき

読み取れる要点は3つです。

  1. 最悪保証が欲しいならKMPかZ-algorithm。どんな入力でも線形が確定します
  2. 実測の速さが欲しいならBoyer-Moore系。理論上の最悪は悪いのに、実務ではもっとも速いことが多いという逆転が起きます
  3. 複数パターンならRabin-Karp(長さが揃う場合)かAho-Corasick(長さがばらばらでも可)

実務ではどうなっているか

ここからが本題かもしれません。現実の高速な文字列検索は、上の教科書的アルゴリズムをそのまま使ってはいません。

Pythonの instr.find は何をしているか

CPythonの部分文字列検索は Objects/stringlib/fastsearch.h に実装されています。ファイル冒頭のコメントには、次のように書かれています。

fast search/count implementation, based on a mix between boyer-
moore and horspool, with a few more bells and whistles on the top.
If the strings are long enough, use Crochemore and Perrin's Two-Way
algorithm, which has worst-case O(n) runtime and best-case O(n/k).

つまり入力の大きさによって使い分けています。ソースの FASTSEARCH 関数を読むと、分岐は次のようになっています。

  1. パターンが1文字なら find_char(内部で memchr に落ちることがある)
  2. n < 2500、または(m < 100 かつ n < 30000)、または m < 6 なら default_find。これがBoyer-Moore / Horspool系にBloomフィルタを足したもので、最悪は O(n * m) です
  3. それ以上の大きさで、パターンがテキストの約33パーセント未満なら、Crochemore-Perrinのtwo-wayアルゴリズムを使う
  4. どちらでもない場合は adaptive_find。素朴に走らせつつ、無駄な部分一致が積み上がってきたらtwo-wayへ切り替える

two-wayアルゴリズムはMaxime CrochemoreとDominique Perrinが1991年に発表したもので("Two-way string-matching", Journal of the ACM 38(3), pp.650-674)、KMPの前向きの照合とBoyer-Mooreの後ろ向きの照合を組み合わせ、前処理 O(m)・検索 O(n)・追加空間 O(log m) を達成します。glibcやmusl、newlibの strstr / memmem でも使われている、実装界の定番です。

Pythonでこれが導入されたのは3.10からで、What's Newには次のように記載されています。

Substring search functions such as str1 in str2 and str2.find(str1)
now sometimes use Crochemore & Perrin's "Two-Way" string searching
algorithm to avoid quadratic behavior on long strings.

貢献者はDennis Sweeney、bpo-41972です。ここで押さえておきたいのは、この最適化は前方検索にのみ適用される点です。ソースを見ると FAST_RSEARCHrfind など)は無条件に default_rfind に落ちており、two-wayは通りません。巨大な文字列に対する rfind は、いまも最悪 O(n * m) になり得ます。

grepとripgrepが使っている手法

GNU grepの作者Mike Haertelは、2010年のメーリングリスト投稿で高速さの理由をこう説明しています。

GNU grep uses the well-known Boyer-Moore algorithm, which looks first
for the final letter of the target string, and uses a lookup table to
tell it how far ahead it can skip in the input whenever it finds a
non-matching character.
GNU grep is fast because it AVOIDS LOOKING AT EVERY INPUT BYTE.

さらに「入力を行に分割しない(大きなバッファのままBoyer-Mooreで走査し、一致した箇所でだけ改行を探す)」「生のシステムコールを使いコピーを避ける」といった、アルゴリズム以外の工夫が並びます。彼のまとめの一文が象徴的です。「プログラムを速くする鍵は、実質的に何もしないようにすることだ」。

複数パターンの場合、GNU grepは kwset.c にCommentz-Walterに近いアルゴリズムを持っています。これはAho-CorasickとBoyer-Mooreを掛け合わせたような手法です。

一方、ripgrepの作者Andrew Gallantは、自身の解説記事でripgrepが取る戦略を説明しています。同記事によれば、Boyer-Mooreで候補位置を見つける部分は memchr のような専用ルーチンに任せるのが有効で、その実装はSIMD命令にコンパイルされて「1回のループで16バイトを調べる」ことができます。ripgrepの工夫は、スキップの起点として末尾の文字ではなく、もっとも出現頻度が低いバイトを選ぶ点にあります。

複数パターンについては、Teddyと呼ばれるSIMDアルゴリズムを使います。同記事の説明では、TeddyはIntelのHyperscanの一部としてGeoffrey Langdaleが考案したもので、16バイトをまとめて比較して候補位置を絞り込みます。リテラルの選択肢が並ぶ正規表現では、正規表現エンジンを使わずAho-Corasickに落とすこともあります。

Rustの memchr クレート(ripgrepと同じ作者)のREADMEはさらに具体的で、memmem の部分文字列検索は小さなテキストにはRabin-Karp、短いパターンには頻度ヒューリスティックを使ったSIMD版、それ以外はtwo-wayアルゴリズムという3段構えだと説明されています。この記事で学んだ手法がそのまま並んでいるのが分かります。

正規表現エンジンとの違い

固定文字列の検索と正規表現マッチングは、目的が近いわりに実装の性質がまったく違います。

  • 固定文字列はスキップできる: パターンに無い文字を踏んだら大きく飛べます。これが劣線形の源です
  • 正規表現は状態機械で1文字ずつ進むのが基本: ただし実用エンジンは「必ず現れるリテラル」を抽出して、まず高速な固定文字列検索で候補を絞り、そこからマッチを試すという最適化を行います

ripgrepの記事によれば、ripgrep、GNU grep、git grepはいずれも有限オートマトンに基づくエンジンを使っており、これは「どんな正規表現・どんな入力でも線形時間で完了する」ことを保証します。一方でバックトラッキング型のエンジンは、この保証がない代わりに後方参照などの機能を持ちます。この違いは正規表現 実践入門でも触れた、いわゆるReDoS(正規表現による処理時間爆発)の背景でもあります。

自分で実装すべきか

結論から言うと、実務でこれらを自作する場面はほとんどありません。判断の目安は次の通りです。

  • 標準ライブラリで済むなら使う: Pythonの in / str.find、Cの strstr / memmem、Javaの String.indexOf などは、上で見た通りすでに高度に最適化されています。素朴に書いたKMPがこれらに勝つことはまずありません
  • ライブラリで足りないケースだけ検討する: 複数パターン同時検索、ストリームを溜めずに照合したい、独自の比較規則(大文字小文字無視、正規化つき)が要る、といった場面です
  • アルゴリズムの知識は「選ぶため」に要る: 自作しないとしても、なぜ rfind だけ遅くなり得るのか、なぜDNA配列の検索が遅いのか、といった判断には理解が不可欠です

複数パターンならAho-Corasick

最後に、単一パターンの枠を超える定番を紹介します。Aho-Corasickアルゴリズムは、Alfred V. AhoとMargaret J. Corasickが1975年に発表しました("Efficient string matching: An aid to bibliographic search", Communications of the ACM 18(6), pp.333-340)。

考え方は「KMPを複数パターンに一般化したもの」と捉えると分かりやすいです。

  1. すべてのパターンを1本のトライ木に格納する
  2. 各ノードにfailureリンク(一致に失敗したとき、そこまでの接尾辞と一致する最長のノードへ飛ぶ辺)を張る
  3. テキストを1文字ずつ流し込み、この有限状態機械を進める

計算量はパターン全体の長さ、テキストの長さ、そして出力される一致の個数の和に対して線形です。パターンが何個あっても、テキストは1回しか走査しません。NGワードフィルタ、ウイルス定義のシグネチャ照合、ログからの複数キーワード抽出といった用途の定番です。歴史的には、Unixの fgrep コマンドがこのアルゴリズムを基礎にしていました。

なお、テキスト側を何度も検索する(パターンが次々変わる)ケースでは、テキスト側を前処理する接尾辞配列接尾辞木が向きます。この記事で扱った手法はすべて「パターンを前処理する」タイプなので、方向性が逆です。用途に応じて、どちら側を前処理するのかを最初に決めるのが設計の分かれ道になります。

まとめ

  • 素朴な全探索は最悪 O(n * m)aaaa...aaaab のように「あと1文字で一致」が繰り返される入力で本当に劣化する
  • Rabin-Karp(Karp and Rabin, 1987年)はローリングハッシュで窓を O(1) 更新し、平均 O(n)。ハッシュ一致は候補にすぎないので必ず実文字列で検証する。長さの揃った複数パターンに強い
  • KMP(Knuth, Morris and Pratt, 1977年)はLPS配列を前処理し、テキストのポインタを一度も戻さないことで前処理 Θ(m) + 検索 Θ(n) を保証する。ストリーム処理に向く
  • Boyer-Moore(Boyer and Moore, 1977年)は後ろから比較し、bad character ruleとgood suffix ruleで大きく飛ぶ。最良 Ω(n / m) の劣線形。Galil ruleなしだと最悪 O(n * m)。簡略版のHorspool(1980年)はbad character ruleのみで実装が短い
  • Z-algorithmはZ配列を O(n) で構築し、パターンと区切り文字とテキストを連結することで検索に使える。区切り文字が両方に現れないことが前提
  • 実務では、CPythonの str.findCrochemore-Perrinのtwo-wayアルゴリズムとHorspool系を入力サイズで切り替え(前方検索のみ、rfind は対象外)、grepはBoyer-Mooreと「1バイトも余計に見ない」設計、ripgrepは memchr のSIMDやTeddyを組み合わせている
  • 複数パターンならAho-Corasick(Aho and Corasick, 1975年)。テキストを1回走査するだけで、全パターンの全出現を線形時間で列挙できる

まずは素朴な全探索とKMPを自分で書き、aaaa...a に対する比較回数を数えてみてください。片方だけが爆発する様子を目の当たりにすると、LPS配列という一見地味な前処理が何を買っているのかが体感できます。そのうえでBoyer-Mooreの「後ろから見る」発想に触れると、理論上の最悪計算量と実測の速さが必ずしも一致しないという、アルゴリズム選択の面白さが見えてくるはずです。

参考リンク

Union-Find(素集合データ構造)入門 - 経路圧縮とunion by sizeで連結性を高速に管理する

Union-Find(素集合データ構造)入門 - 経路圧縮とunion by sizeで連結性を高速に管理する

32

要素のグループ分けを管理するUnion-Find(素集合データ構造 / DSU)を解説します。find・union・sameという3つの基本操作、森による表現、素朴な実装がO(n)に劣化する理由、経路圧縮とunion by sizeという2つの最適化、Tarjanが示したならしO(α(n))という計算量、クラスカル法・連結成分の数え上げ・グリッドの島の数・サイクル検出・重み付きUnion-Findといった応用、各言語の標準ライブラリ事情まで、実際に動くPythonコードで整理します。

ヒープと優先度付きキュー入門 - O(log n)で最小・最大を取り出すデータ構造

ヒープと優先度付きキュー入門 - O(log n)で最小・最大を取り出すデータ構造

18

優先度付きキューという抽象データ型と、その代表的な実装である二分ヒープを解説。完全二分木の配列表現、sift-up/sift-down、push/popがO(log n)・build-heapがO(n)といった計算量、ヒープソート、各言語の標準ライブラリ、Dijkstraやハフマン符号化・Top-K・中央値の2ヒープ法まで、Pythonコードで整理します。