Bit-Parallel Approximate Matching with Packed DP Rows

Approximate matching asks where a pattern occurs in a text allowing up to k edits — insertions, deletions, substitutions. The standard dynamic programme fills an (m+1)×(n+1) table of edit distances between every prefix of the pattern and every prefix of the text, which costs O(nm) time and, if only the final row is kept, O(m) space. For a pattern of 64 characters that is 64 cells per text character, each cell an add, two compares, and a minimum.

The bit-parallel formulation replaces the 64 cells with 64 bits of a single register and computes the whole column with a handful of bitwise operations. It is not an approximation of the dynamic programme; it computes the same values, encoded differently. The encoding replaces m cell updates per text character with a fixed number of word-sized operations, so the work per character no longer grows with the pattern length.

Encoding a column as bits

The dynamic programme’s column j holds, for each i, the edit distance D[i][j] between the pattern prefix of length i and the text prefix of length j. Two properties of adjacent cells make the encoding possible. The values change by at most one between horizontally or vertically adjacent cells, so the difference D[i][j] − D[i−1][j] is −1, 0, or 1 in the standard formulation: adding one more character to the pattern prefix can raise the distance to a fixed text prefix by one, leave it unchanged, or lower it by one when the new character lets an existing alignment improve.

So a column can be represented by a base value and a bit vector of the vertical differences. Myers’ formulation goes further and encodes, for each row, whether the horizontal difference is +1, 0, or −1, using two bit vectors: one marking cells where the value increases relative to the cell to the left, one where it decreases. The current column’s distances are then recovered by a prefix sum of ±1 steps, and the score for a match ending at the current text position is the value at row m, maintained as an integer counter rather than recomputed.

The update uses three bit vectors, conventionally named for their role:

• Pv, positive vertical differences: rows where D[i][j] − D[i−1][j] = 1.
• Mv, negative vertical differences: rows where D[i][j] − D[i−1][j] = −1, the pattern prefix having moved one edit closer to the text prefix than the shorter prefix above it.
• Eq, the match mask for the current text character: bit i is set if pattern[i] equals the character.

The recurrence computes the new Pv, the new Mv, and updates the score, all with about fifteen bitwise operations per text character. The bit at position m−1 of the intermediate Ph and Mh masks says whether the score increases or decreases, which is how the running distance is kept without touching the whole column.

#include <stdint.h>
#include <string.h>

/* Match mask for each byte value: bit i set iff pattern[i] == c. */
static void build_peq(const char *pat, size_t m, uint64_t peq[256]) {
    memset(peq, 0, 256 * sizeof(uint64_t));
    for (size_t i = 0; i < m; i++)
        peq[(unsigned char)pat[i]] |= (uint64_t)1 << i;
}

/* Semi-global: pattern matched against any substring of text.
   Returns the smallest edit distance found, and its end offset via *end. */
static unsigned search(const char *pat, size_t m, const char *text, size_t n,
                       uint64_t *end) {
    uint64_t peq[256];
    build_peq(pat, m, peq);

    const uint64_t high = (uint64_t)1 << (m - 1);
    uint64_t pv = ~(uint64_t)0;   /* all ones: column 0 of the DP table */
    uint64_t mv = 0;
    unsigned score = (unsigned)m, best = (unsigned)m;

    for (size_t j = 0; j < n; j++) {
        uint64_t eq = peq[(unsigned char)text[j]];
        uint64_t xv = eq | mv;
        uint64_t xh = (((eq & pv) + pv) ^ pv) | eq;
        uint64_t ph = mv | ~(xh | pv);
        uint64_t mh = pv & xh;

        if (ph & high) score++;
        else if (mh & high) score--;

        if (score < best) { best = score; if (end) *end = j; }

        ph = ph << 1;             /* no seed: row 0 has zero horizontal delta here */
        mh = mh << 1;
        pv = mh | ~(xv | ph);
        mv = ph & xv;
    }
    return best;
}

The initial pv = ~0 is not arbitrary. Column 0 of the table is D[i][0] = i, so every vertical difference is +1, and Pv starts with all bits set. The ph << 1 shift moves the horizontal deltas one row up, and the missing low bit is the horizontal delta out of row 0, which is fixed by the boundary condition: zero in the semi-global variant, where an alignment may start anywhere, and one in the global variant, where the distance from an empty pattern to a text prefix of length j is j. Writing (ph << 1) | 1 in the semi-global case feeds the value for row 0 in the global variant, charging one edit per text position before any pattern character has been read and inflating the score; the bug is subtle because the algorithm still runs and still reports plausible numbers.

The ph & high test reads the top bit of the word, which tracks row m, the full pattern. A set bit there means the distance to the pattern prefix increased at this text position; a set bit in mh means it decreased. That is the entire mechanism for tracking the score, and it is why the algorithm needs no table.

What semi-global means here

The code above starts with a full row of ones in Pv, which is the encoding of D[i][0] = i, and suppresses the seed bit so that row 0 stays at zero. Row 0 free means the pattern may begin at any text position without cost, so the score reports the best alignment ending anywhere rather than only one ending at the last text character. That is the variant wanted for searching. The global edit distance is the other initialisation: D[0][j] = j is expressed by seeding the low bit with one on every step, and the answer is then the single cell at the bottom-right of the table. The difference is a one-bit constant per step, not a change in the recurrence itself.

For a search with a maximum error budget k, the algorithm can stop reporting once the score exceeds k, and the score is monotone in the sense that it changes by at most one per character. Early termination on a long text is therefore possible, but it requires the pattern to have failed; the score can also drop again, so termination must be based on the best score over a window rather than the current value.

The same bit-vector machinery computes exact matching as the special case k = 0: the score starts at m and never reaches 0 unless a full match occurred, and the Shift-And algorithm is exactly this recurrence with the score tracking removed.

Shift-And as the degenerate case

For exact matching, the state is one bit vector D where bit i is set if the pattern prefix of length i+1 matches a suffix of the text read so far. The update is

D = ((D << 1) | 1) & Eq(c)

A match is reported when bit m−1 of D is set. The | 1 seeds the first row, which always matches the empty prefix; the shift propagates each existing match one character to the right; the AND with the match mask kills every position where the character does not agree.

This is a simulation of the nondeterministic automaton of all pattern prefixes, run in parallel across the bits of a word. Its cost is one shift, one OR, one AND, and one table lookup per text character, independent of the pattern length up to the word size. For a pattern of up to 64 bytes on a 64-bit machine it examines each text character with a constant number of operations, and the approximate version generalises the same recurrence.

Shift-Or is the same algorithm with the match mask complemented, so a set bit means mismatch rather than match, and the accepting condition inverted: a zero bit at the top signals a match. The two names are used interchangeably.

The word-size wall and block decomposition

The bit-parallel recurrence processes m rows, so it works unmodified only when m ≤ w, the machine word size. Two strategies extend it.

Block decomposition splits the pattern into ⌈m/w⌉ blocks and runs one bit-vector state per block, carrying a match from one block to the next. For exact matching this is straightforward: block b’s state is advanced by the same character, and a match in block b at bit w−1 seeds bit 0 of block b+1. The carry costs an extra operation per block per character, so the total cost is ⌈m/w⌉ operations per text character, which is the same asymptotic complexity as the table method but with a word-sized constant. The catch is that the carry must be exact: a match spanning a block boundary must be detected at the moment the last character of the boundary-crossing prefix arrives, and the ordering of seeds between blocks must respect the direction of propagation.

The approximate case blocks differently. The score is a single integer for the whole pattern, and block decomposition requires each block’s local score plus an adjustment for the boundary, which is not exact in general. A blocked implementation therefore treats the per-block result as a bound that can exceed the true distance, and uses it for filtering: a block that reports a distance above k prunes the alignment, and surviving candidates are verified. This is the standard shape of approximate matchers — a fast probabilistic or bounded filter followed by exact verification of a small candidate set.

Pattern length	Strategy	Work per text character
≤ 64	single word, Myers	~15 bitwise ops
≤ 128	two words with carry	~30 ops
m ≫ 64	blocked filter plus verification	⌈m/64⌉ ops, plus verification

Why the bit-vector form is faster in practice

The table method performs a minimum of three values per cell, with data-dependent branches that mispredict. The bit-parallel method performs a fixed sequence of register operations with no branches except the score update, and the whole state lives in two registers. On a modern superscalar core the operations issue in parallel and the dependency chain is short, so the per-character cost is a handful of cycles rather than a memory round trip per cell.

The advantage disappears when the alphabet is large and the match mask does not fit in cache. The peq table has one word per alphabet symbol: 2 KB for a byte alphabet, which fits in L1, but 512 KB for a 16-bit alphabet and more than any cache for code points. For Unicode patterns the table must be built per pattern and touched sparsely, and the lookup becomes the bottleneck. Byte-oriented matching, or a map from character to a small index, restores the locality.

A second limitation is that the recurrence computes edit distance, not arbitrary scoring. Affine gap costs, where opening a gap is more expensive than extending one, need additional bit vectors to track gap state, and the constant grows. Position-specific scoring matrices, as used in biological alignment, do not fit the bit encoding at all because the cost of a substitution depends on which characters are involved, not merely whether they match.

The third limitation is k. The bit-vector formulation computes exact distances, so the error budget is a reporting threshold rather than a bound on work. Algorithms that exploit a small k to prune — the classic Wu–Manber filtering, or the pigeonhole principle that splits a pattern with k errors into k+1 pieces each of which must match exactly — reduce work further when k is small and m is large, at the cost of a verification phase. The bit-parallel core is the verification primitive those filters call.

Where the other matchers fit

Exact matching by bit-parallel simulation is the same automaton as the prefix-function method, with the state held as bits instead of an integer. The KMP state says how long the current match is; the Shift-And state says, for every prefix length, whether it matches, which is more information packed into fewer operations. Neither extends to a dictionary of patterns, which is Aho–Corasick territory, though a bit-parallel simulation of a trie’s NFA is possible and is what some regular-expression engines do for small patterns.

Approximate matching over a large text with many queries is better served by an index. The suffix automaton can be run with a dynamic programme over its states and the pattern in O(states × |P|) time, and the FM-index supports backward search with a bounded number of mismatches by branching on the alphabet at each step. Those structures pay a construction cost proportional to the text; the bit-parallel matcher pays nothing and scans linearly, which is why it remains the right choice for a single pass over a large text with a modest pattern.