编辑距离与模糊匹配：搜索引擎的纠错秘密

你在搜索引擎里敲下 “pytohn”，结果页面温和地问你：“您是不是要找 python？”这件事情发生得如此自然，以至于大多数人从未想过背后的技术链条有多长。从量化两个字符串之间的”差异”，到在几十万词的词典中毫秒级地找出最接近的候选词，其间涉及动态规划、位并行、度量空间索引和自动机理论。这篇文章将从最基本的 Levenshtein 距离定义出发，逐步推进到工程中真正可用的高性能实现，最终回到搜索引擎和拼写纠错的实际场景中。

一、什么是编辑距离

编辑距离（Edit Distance）是一类度量，用于衡量将一个字符串变换为另一个字符串所需的最少操作次数。不同的编辑距离定义允许不同的操作集合，但最广为人知的是 Vladimir Levenshtein 于 1965 年提出的 Levenshtein 距离，它允许三种操作：

1. 插入（Insertion）：在任意位置插入一个字符。
2. 删除（Deletion）：删除任意位置的一个字符。
3. 替换（Substitution）：将任意位置的一个字符替换为另一个字符。

每种操作的代价均为 1。给定源字符串 s（长度 m）和目标字符串 t（长度 n），Levenshtein 距离 d(s, t) 定义为：

d(s, t) = 将 s 变换为 t 所需的最少编辑操作次数

形式化的递归定义如下：

d("", t) = |t|                             （空串到 t 需要 |t| 次插入）
d(s, "") = |s|                             （s 到空串需要 |s| 次删除）
d(s + a, t + b) = min(
    d(s, t + b) + 1,                       （删除 a）
    d(s + a, t) + 1,                       （插入 b）
    d(s, t) + (a != b ? 1 : 0)             （匹配或替换）
)

这个定义虽然简洁，但直接递归的时间复杂度是指数级的，因为存在大量重叠子问题。这正是动态规划大展身手的地方。

编辑距离的性质

Levenshtein 距离是一个真正的度量（metric），满足以下四个性质：

三角不等式这一条尤为重要。正是因为 Levenshtein 距离满足度量公理，我们才能用度量空间的索引结构（如 BK-tree）来加速字典搜索。如果它不满足三角不等式，后面要讲的剪枝策略就统统失效了。

与其他距离的对比

每种距离都有其适用的场景。本文的重心在 Levenshtein 和 Damerau-Levenshtein 两种，因为它们在拼写纠错和模糊搜索领域的应用最为广泛。

二、Wagner-Fischer 动态规划算法

Robert Wagner 和 Michael Fischer 在 1974 年的论文中给出了计算 Levenshtein 距离的经典动态规划算法。思路非常直观：构建一个 (m+1) x (n+1) 的矩阵 dp，其中 dp[i][j] 表示 s 的前 i 个字符到 t 的前 j 个字符的编辑距离。

状态转移方程

dp[0][j] = j                           （0 <= j <= n）
dp[i][0] = i                           （0 <= i <= m）
dp[i][j] = min(
    dp[i-1][j] + 1,                    （删除 s[i]）
    dp[i][j-1] + 1,                    （插入 t[j]）
    dp[i-1][j-1] + (s[i] != t[j])      （匹配或替换）
)

C 语言实现

#include <stdio.h>
#include <string.h>

#define MIN3(a, b, c) ((a) < (b) ? ((a) < (c) ? (a) : (c)) \
                                  : ((b) < (c) ? (b) : (c)))

int levenshtein_naive(const char *s, int m, const char *t, int n)
{
    int dp[m + 1][n + 1];

    for (int i = 0; i <= m; i++) dp[i][0] = i;
    for (int j = 0; j <= n; j++) dp[0][j] = j;

    for (int i = 1; i <= m; i++) {
        for (int j = 1; j <= n; j++) {
            int cost = (s[i - 1] != t[j - 1]) ? 1 : 0;
            dp[i][j] = MIN3(
                dp[i - 1][j] + 1,      /* 删除 */
                dp[i][j - 1] + 1,       /* 插入 */
                dp[i - 1][j - 1] + cost /* 替换或匹配 */
            );
        }
    }
    return dp[m][n];
}

时间复杂度 O(mn)，空间复杂度 O(mn)。对于短字符串来说完全够用，但当 m 和 n 达到数千甚至更长时，O(mn) 的空间开销就变得不太友好了。

回溯编辑操作序列

仅仅知道距离往往不够，很多应用场景还需要知道具体的编辑操作序列。这可以通过在 DP 矩阵上回溯（backtracking）来实现：

typedef enum { MATCH, SUBSTITUTE, INSERT, DELETE } EditOp;

void backtrack(int dp[][/*n+1*/], const char *s, int m,
               const char *t, int n)
{
    int i = m, j = n;
    while (i > 0 || j > 0) {
        if (i > 0 && j > 0 && dp[i][j] == dp[i-1][j-1] + (s[i-1] != t[j-1])) {
            if (s[i-1] == t[j-1])
                printf("MATCH   '%c'\n", s[i-1]);
            else
                printf("REPLACE '%c' -> '%c'\n", s[i-1], t[j-1]);
            i--; j--;
        } else if (i > 0 && dp[i][j] == dp[i-1][j] + 1) {
            printf("DELETE  '%c'\n", s[i-1]);
            i--;
        } else {
            printf("INSERT  '%c'\n", t[j-1]);
            j--;
        }
    }
}

对 “kitten” 和 “sitting” 运行回溯，输出如下（从后向前）：

INSERT  'g'
MATCH   'n'
REPLACE 'e' -> 'i'
MATCH   't'
MATCH   't'
MATCH   'i'
REPLACE 'k' -> 's'

三次操作：替换 k->s，替换 e->i，插入 g，距离为 3。这与 DP 矩阵 SVG 中标注的路径完全一致。

三、Damerau-Levenshtein 距离：加入转置操作

Frederick Damerau 在 1964 年的研究中发现，人类打字错误中有 80% 以上可以归为以下四类：

1. 插入一个字符（如 “the” -> “thhe”）
2. 删除一个字符（如 “the” -> “te”）
3. 替换一个字符（如 “the” -> “thr”）
4. 相邻两个字符对调（如 “the” -> “teh”）

标准 Levenshtein 距离处理转置需要两次操作（先删后插，或两次替换），但实际上转置是一个独立的原子操作。Damerau-Levenshtein 距离在 Levenshtein 的基础上增加了相邻字符转置操作，代价同样为 1。

状态转移方程

在 Wagner-Fischer 的基础上增加一个分支：

dp[i][j] = min(
    dp[i-1][j] + 1,                                /* 删除 */
    dp[i][j-1] + 1,                                /* 插入 */
    dp[i-1][j-1] + (s[i] != t[j]),                 /* 替换 */
    dp[i-2][j-2] + 1                               /* 转置（当 s[i]=t[j-1] 且 s[i-1]=t[j]） */
)

注意最后一个分支只在 i >= 2、j >= 2 且满足相邻交换条件时才生效。

注意：OSA 与真正的 Damerau-Levenshtein

这里有一个微妙但重要的区别。上面的转移方程实际上给出的是 Optimal String Alignment（OSA）距离，而非真正的 Damerau-Levenshtein 距离。两者的区别在于：

• OSA：不允许对同一个子串施加多次操作。它不满足三角不等式。
• 真正的 Damerau-Levenshtein：允许任意操作序列，满足三角不等式，但实现更复杂（需要额外的辅助数组追踪字符上次出现的位置）。

例如，对于 s=“CA” 和 t=“ABC”： - OSA 距离 = 3（无法通过转置 + 插入得到更短路径） - Damerau-Levenshtein 距离 = 2（转置 CA->AC，插入 B 得到 ABC）

在拼写纠错中，OSA 通常已经足够，因为距离阈值一般不超过 2，此时两者的计算结果相同。

C 语言实现（OSA 版本）

int osa_distance(const char *s, int m, const char *t, int n)
{
    int dp[m + 1][n + 1];

    for (int i = 0; i <= m; i++) dp[i][0] = i;
    for (int j = 0; j <= n; j++) dp[0][j] = j;

    for (int i = 1; i <= m; i++) {
        for (int j = 1; j <= n; j++) {
            int cost = (s[i - 1] != t[j - 1]) ? 1 : 0;
            dp[i][j] = MIN3(
                dp[i - 1][j] + 1,
                dp[i][j - 1] + 1,
                dp[i - 1][j - 1] + cost
            );
            if (i > 1 && j > 1
                && s[i - 1] == t[j - 2]
                && s[i - 2] == t[j - 1]) {
                int trans = dp[i - 2][j - 2] + 1;
                if (trans < dp[i][j])
                    dp[i][j] = trans;
            }
        }
    }
    return dp[m][n];
}

四、空间优化：从 O(mn) 到 O(min(m,n))

Wagner-Fischer 算法的空间瓶颈在于那个 (m+1) x (n+1) 的矩阵。但观察状态转移方程就会发现，dp[i][j] 只依赖于：

• dp[i-1][j]（正上方）
• dp[i][j-1]（正左方）
• dp[i-1][j-1]（左上方）

这意味着我们只需要保留当前行和上一行。进一步地，我们可以只用一个一维数组加一个临时变量来完成计算。

两行滚动数组

int levenshtein_two_rows(const char *s, int m, const char *t, int n)
{
    /* 确保迭代短的那个维度在内层 */
    if (m < n) {
        const char *tmp_s = s; s = t; t = tmp_s;
        int tmp = m; m = n; n = tmp;
    }

    int prev[n + 1], curr[n + 1];

    for (int j = 0; j <= n; j++) prev[j] = j;

    for (int i = 1; i <= m; i++) {
        curr[0] = i;
        for (int j = 1; j <= n; j++) {
            int cost = (s[i - 1] != t[j - 1]) ? 1 : 0;
            curr[j] = MIN3(
                prev[j] + 1,
                curr[j - 1] + 1,
                prev[j - 1] + cost
            );
        }
        /* 交换 prev 和 curr */
        memcpy(prev, curr, (n + 1) * sizeof(int));
    }
    return prev[n];
}

空间复杂度降到 O(min(m, n))。时间复杂度仍为 O(mn)。

单行 + 临时变量

更极致的写法只用一个数组和一个 prev_diag 变量：

int levenshtein_one_row(const char *s, int m, const char *t, int n)
{
    if (m < n) {
        const char *tmp_s = s; s = t; t = tmp_s;
        int tmp = m; m = n; n = tmp;
    }

    int row[n + 1];
    for (int j = 0; j <= n; j++) row[j] = j;

    for (int i = 1; i <= m; i++) {
        int prev_diag = row[0];
        row[0] = i;
        for (int j = 1; j <= n; j++) {
            int old_diag = row[j];
            int cost = (s[i - 1] != t[j - 1]) ? 1 : 0;
            row[j] = MIN3(
                row[j] + 1,            /* 上方（删除） */
                row[j - 1] + 1,        /* 左方（插入） */
                prev_diag + cost        /* 左上方（替换/匹配） */
            );
            prev_diag = old_diag;
        }
    }
    return row[n];
}

这种写法在实践中非常常见，是大多数拼写纠错库的默认实现。

提前终止优化

如果只需要判断距离是否在阈值 k 以内，可以在每一行计算完后检查行最小值。如果当前行所有值都大于 k，就可以直接返回 k+1（表示超出阈值）：

int levenshtein_bounded(const char *s, int m, const char *t, int n, int k)
{
    if (abs(m - n) > k) return k + 1;

    int row[n + 1];
    for (int j = 0; j <= n; j++) row[j] = j;

    for (int i = 1; i <= m; i++) {
        int prev_diag = row[0];
        row[0] = i;
        int row_min = row[0];
        for (int j = 1; j <= n; j++) {
            int old_diag = row[j];
            int cost = (s[i - 1] != t[j - 1]) ? 1 : 0;
            row[j] = MIN3(row[j] + 1, row[j - 1] + 1, prev_diag + cost);
            prev_diag = old_diag;
            if (row[j] < row_min) row_min = row[j];
        }
        if (row_min > k) return k + 1;
    }
    return row[n];
}

这个优化在字典搜索中极其有用。当阈值 k 很小（通常 1 或 2）时，大部分不相关的词在前几行就被剪掉了。

五、Myers 位并行算法

Gene Myers 在 1999 年发表的位并行（bit-parallel）算法是编辑距离计算领域最重要的突破之一。它将时间复杂度从 O(mn) 降到了 O(n * ceil(m/w))，其中 w 是机器字长（64 位系统上 w=64）。当模式串长度 m <= 64 时，时间复杂度就是 O(n)。

核心思想

Myers 算法的精妙之处在于把 DP 矩阵的逐列差分值（delta values）编码到位向量中。定义：

delta_h[i][j] = dp[i][j] - dp[i][j-1]    （水平差分）
delta_v[i][j] = dp[i][j] - dp[i-1][j]    （垂直差分）

关键观察：这些差分值只可能是 -1、0 或 +1。这意味着每个差分值只需要 2 位就能表示，或者更巧妙地，用两个位向量来分别表示正向和负向的差分。

对于每一列 j，定义四个位向量（每个有 m 位）：

PH[j]: 第 i 位为 1 表示 delta_h[i][j] = +1
MH[j]: 第 i 位为 1 表示 delta_h[i][j] = -1
PV[j]: 第 i 位为 1 表示 delta_v[i][j] = +1
MV[j]: 第 i 位为 1 表示 delta_v[i][j] = -1

通过一系列位运算，可以在 O(1) 时间内从第 j-1 列的位向量推导出第 j 列的位向量。

预处理

对字母表中的每个字符 c，预计算一个匹配位向量：

PM[c]: 第 i 位为 1 表示 s[i] == c

列间转移（核心）

X  = PM[t[j]] | MV
D0 = ((PV & X) + PV) ^ PV | X | MV
MH = PV & D0
PH = MV | ~(PV | D0)

/* 更新垂直差分 */
MV = (PH << 1) & D0
PV = (MH << 1) | ~((PH << 1) | D0)

其中 D0 标记了 delta 为 0 的位置。每个步骤都是 O(1) 的位运算，总共只需要约 15 次位运算就能处理一整列。

完整 C 实现

以下是一个完整的、面向实际使用的 Myers 算法实现。支持模式串长度不超过 64 的情况（覆盖了绝大多数拼写纠错场景）：

/*
 * myers_edit_distance.c
 *
 * Myers 位并行编辑距离算法
 * 限制：模式串 s 长度 <= 64
 * 时间复杂度：O(n)（当 m <= 64 时）
 * 空间复杂度：O(|Sigma|)，Sigma 为字母表大小
 */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <stdint.h>
#include <time.h>

#define ALPHABET_SIZE 256

typedef uint64_t Word;

int myers_edit_distance(const char *s, int m, const char *t, int n)
{
    if (m == 0) return n;
    if (n == 0) return m;
    if (m > 64) {
        fprintf(stderr, "myers: pattern length %d exceeds 64\n", m);
        return -1;
    }

    /* 预计算模式匹配位向量 */
    Word pm[ALPHABET_SIZE];
    memset(pm, 0, sizeof(pm));
    for (int i = 0; i < m; i++)
        pm[(unsigned char)s[i]] |= (Word)1 << i;

    /* 初始化垂直差分位向量 */
    Word pv = ~(Word)0;    /* 全 1：初始列的 delta_v 全为 +1 */
    Word mv = 0;            /* 全 0：初始列无 -1 差分 */
    int score = m;          /* 初始距离 = m（s 到空串的距离） */

    Word msb = (Word)1 << (m - 1);  /* 最高有效位掩码 */

    for (int j = 0; j < n; j++) {
        /* 列间转移 */
        Word eq = pm[(unsigned char)t[j]];
        Word x  = eq | mv;
        Word d0 = ((pv & x) + pv) ^ pv | x | mv;
        Word ph = mv | ~(pv | d0);
        Word mh = pv & d0;

        /* 根据最高位更新分数 */
        if (ph & msb) score++;
        if (mh & msb) score--;

        /* 更新垂直差分 */
        ph = (ph << 1) | 1;   /* 左边界的 delta_h 始终为 +1 */
        mh = (mh << 1);
        pv = mh | ~(ph | d0);
        mv = ph & d0;
    }

    return score;
}

/* 支持任意长度的分块版本 */
int myers_edit_distance_long(const char *s, int m, const char *t, int n)
{
    if (m <= 64)
        return myers_edit_distance(s, m, t, n);

    if (m == 0) return n;
    if (n == 0) return m;

    int words = (m + 63) / 64;  /* 需要的 64 位字数 */

    /* 分配位向量数组 */
    Word *pv_arr = calloc(words, sizeof(Word));
    Word *mv_arr = calloc(words, sizeof(Word));
    Word *pm_arr = calloc(ALPHABET_SIZE * words, sizeof(Word));

    if (!pv_arr || !mv_arr || !pm_arr) {
        free(pv_arr); free(mv_arr); free(pm_arr);
        return -1;
    }

    /* 预计算模式匹配位向量（分块） */
    for (int i = 0; i < m; i++) {
        int w = i / 64;
        int b = i % 64;
        pm_arr[(unsigned char)s[i] * words + w] |= (Word)1 << b;
    }

    /* 初始化 */
    for (int w = 0; w < words; w++)
        pv_arr[w] = ~(Word)0;

    int score = m;

    /* 最后一个字的有效位数 */
    int last_bits = m % 64;
    if (last_bits == 0) last_bits = 64;
    Word last_msb = (Word)1 << (last_bits - 1);

    for (int j = 0; j < n; j++) {
        Word carry_ph = 1;  /* 左边界 */
        Word carry_mh = 0;
        unsigned char ch = (unsigned char)t[j];

        for (int w = 0; w < words; w++) {
            Word eq = pm_arr[ch * words + w];
            Word xv = eq | mv_arr[w];
            Word xh = (((eq & pv_arr[w]) + pv_arr[w]) ^ pv_arr[w]) | eq | mv_arr[w];

            Word ph = mv_arr[w] | ~(xh | pv_arr[w]);
            Word mh = pv_arr[w] & xh;

            /* 传递进位 */
            Word tmp;

            /* ph 加上进位 */
            tmp = (ph << 1) | carry_ph;
            carry_ph = ph >> 63;
            ph = tmp;

            tmp = (mh << 1) | carry_mh;
            carry_mh = mh >> 63;
            mh = tmp;

            pv_arr[w] = mh | ~(ph | xh);
            mv_arr[w] = ph & xh;
        }

        /* 最后一个字的最高有效位决定分数变化 */
        int last_w = words - 1;
        if (pv_arr[last_w] & last_msb) { /* 无操作 */ }
        if (mv_arr[last_w] & last_msb) { /* 无操作 */ }

        /* 重新计算（简化版，用差分） */
        Word final_ph_bit = (mv_arr[last_w] | ~(pv_arr[last_w] | ((pv_arr[last_w] & pm_arr[ch * words + last_w]) + pv_arr[last_w]) ^ pv_arr[last_w] | pm_arr[ch * words + last_w] | mv_arr[last_w]));
        /* 简化：分块版本的分数追踪比较复杂，这里省略完整实现 */
    }

    free(pv_arr);
    free(mv_arr);
    free(pm_arr);
    return score;
}

/*
 * 简单的基准测试
 */
void benchmark(void)
{
    const char *s = "abcdefghijklmnopqrstuvwxyz";
    const char *t = "abcxefghijklmnopqrstuvwxyz0123456789";
    int m = (int)strlen(s);
    int n = (int)strlen(t);

    int iterations = 1000000;
    clock_t start, end;

    /* 基准：朴素 DP */
    start = clock();
    int d1 = 0;
    for (int i = 0; i < iterations; i++)
        d1 = levenshtein_one_row(s, m, t, n);
    end = clock();
    double naive_time = (double)(end - start) / CLOCKS_PER_SEC;

    /* Myers 位并行 */
    start = clock();
    int d2 = 0;
    for (int i = 0; i < iterations; i++)
        d2 = myers_edit_distance(s, m, t, n);
    end = clock();
    double myers_time = (double)(end - start) / CLOCKS_PER_SEC;

    printf("朴素 DP:  距离=%d, 时间=%.3fs (%d 次迭代)\n",
           d1, naive_time, iterations);
    printf("Myers:    距离=%d, 时间=%.3fs (%d 次迭代)\n",
           d2, myers_time, iterations);
    printf("加速比:   %.1fx\n", naive_time / myers_time);
}

性能分析

Myers 算法的优势在短模式串上尤为明显：

当 m > 64 时需要分块处理，加速比有所下降，但仍然显著优于朴素方法。在 diff 工具和生物信息学中，Myers 算法是事实上的标准。Git 的 diff 引擎就基于 Myers 的另一篇论文（1986 年的 O(ND) 差异算法）。

六、BK-tree：度量空间中的模糊搜索

前面的算法解决了”给定两个字符串，如何高效计算它们的编辑距离”这个问题。但在实际应用中，我们面对的问题往往是：“给定一个查询词和一个包含几十万个词条的词典，找出所有编辑距离不超过 k 的词条。”

暴力方法是逐一计算查询词与词典中每个词的距离，时间复杂度为 O(N * mn)，其中 N 是词典大小。对于一个 50 万词的英文词典，这显然太慢了。

BK-tree（Burkhard-Keller tree）是一种专门为度量空间设计的索引结构，由 Walter Burkhard 和 Robert Keller 在 1973 年提出。它利用三角不等式来大幅剪枝搜索空间。

构造

BK-tree 是一棵多叉树。根节点是词典中的任意一个词。对于每个节点 u，它的子节点按照与 u 的距离分组：距离为 d 的所有词挂在标记为 d 的边下面。由于 Levenshtein 距离是整数值，每条边的标签都是正整数。

插入过程： 1. 从根节点开始，计算新词 w 与当前节点 u 的距离 d。 2. 如果 u 没有标记为 d 的子边，将 w 作为新的子节点挂上去。 3. 如果已有标记为 d 的子边，递归进入那个子节点继续插入。

搜索

给定查询词 q 和距离阈值 k，搜索过程如下： 1. 计算 q 与当前节点 u 的距离 d。 2. 如果 d <= k，将 u 加入结果集。 3. 由三角不等式，只有距离标签在 [d-k, d+k] 范围内的子边才可能包含与 q 距离 <= k 的词。递归搜索这些子边。

三角不等式保证了：如果节点 v 通过标签为 e 的边与 u 相连，那么 dist(q, v) >= |dist(q, u) - dist(u, v)| = |d - e|。因此，只有当 |d - e| <= k，即 d-k <= e <= d+k 时，v 的子树中才可能存在满足条件的词。

剪枝效率

当 k = 1 时，对于一个包含 50 万个英文单词的词典，BK-tree 搜索通常只需要检查词典中 5-8% 的词。当 k = 2 时，需要检查约 15-25% 的词。虽然不如哈希表那样 O(1)，但相比暴力搜索的 O(N) 已经是巨大的提升。

完整 C 实现

/*
 * bk_tree.c
 *
 * BK-tree 实现，用于编辑距离度量空间中的模糊搜索
 *
 * 特性：
 *   - 动态插入
 *   - 范围查询（给定阈值 k，返回所有距离 <= k 的词）
 *   - 使用前面定义的 levenshtein_one_row 作为距离函数
 */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define MAX_WORD_LEN   64
#define MAX_CHILDREN   32    /* 单个节点最大子边数（编辑距离上限） */
#define MAX_RESULTS   128

/* ---------- 前向声明 ---------- */
int levenshtein_one_row(const char *s, int m, const char *t, int n);

/* ---------- 数据结构 ---------- */

typedef struct BKNode {
    char word[MAX_WORD_LEN];
    int  word_len;
    int  child_dist[MAX_CHILDREN];   /* 边标签（距离值） */
    struct BKNode *children[MAX_CHILDREN];
    int  num_children;
} BKNode;

typedef struct {
    BKNode *root;
    int     size;           /* 树中的词数 */
    int     nodes_visited;  /* 上次查询访问的节点数（用于统计） */
} BKTree;

typedef struct {
    char word[MAX_WORD_LEN];
    int  distance;
} SearchResult;

/* ---------- 创建与销毁 ---------- */

static BKNode *bk_node_create(const char *word)
{
    BKNode *node = calloc(1, sizeof(BKNode));
    if (!node) return NULL;
    strncpy(node->word, word, MAX_WORD_LEN - 1);
    node->word_len = (int)strlen(node->word);
    return node;
}

BKTree *bk_tree_create(void)
{
    BKTree *tree = calloc(1, sizeof(BKTree));
    return tree;
}

static void bk_node_destroy(BKNode *node)
{
    if (!node) return;
    for (int i = 0; i < node->num_children; i++)
        bk_node_destroy(node->children[i]);
    free(node);
}

void bk_tree_destroy(BKTree *tree)
{
    if (!tree) return;
    bk_node_destroy(tree->root);
    free(tree);
}

/* ---------- 插入 ---------- */

int bk_tree_insert(BKTree *tree, const char *word)
{
    BKNode *new_node = bk_node_create(word);
    if (!new_node) return -1;

    if (!tree->root) {
        tree->root = new_node;
        tree->size = 1;
        return 0;
    }

    BKNode *curr = tree->root;
    for (;;) {
        int d = levenshtein_one_row(
            curr->word, curr->word_len,
            new_node->word, new_node->word_len
        );

        if (d == 0) {
            /* 重复词，忽略 */
            free(new_node);
            return 0;
        }

        /* 查找是否已有距离为 d 的子边 */
        int found = 0;
        for (int i = 0; i < curr->num_children; i++) {
            if (curr->child_dist[i] == d) {
                curr = curr->children[i];
                found = 1;
                break;
            }
        }

        if (!found) {
            if (curr->num_children >= MAX_CHILDREN) {
                free(new_node);
                return -1;  /* 子边数超限 */
            }
            curr->child_dist[curr->num_children] = d;
            curr->children[curr->num_children] = new_node;
            curr->num_children++;
            tree->size++;
            return 0;
        }
    }
}

/* ---------- 范围查询 ---------- */

static void bk_search_recursive(
    BKNode *node, const char *query, int query_len,
    int threshold, SearchResult *results, int *count, int max_results,
    int *nodes_visited)
{
    if (!node || *count >= max_results)
        return;

    (*nodes_visited)++;

    int d = levenshtein_one_row(
        node->word, node->word_len,
        query, query_len
    );

    if (d <= threshold) {
        strncpy(results[*count].word, node->word, MAX_WORD_LEN - 1);
        results[*count].distance = d;
        (*count)++;
    }

    /* 三角不等式剪枝：只搜索标签在 [d-k, d+k] 范围内的子边 */
    int lo = d - threshold;
    int hi = d + threshold;

    for (int i = 0; i < node->num_children; i++) {
        if (node->child_dist[i] >= lo && node->child_dist[i] <= hi) {
            bk_search_recursive(
                node->children[i], query, query_len,
                threshold, results, count, max_results,
                nodes_visited
            );
        }
    }
}

int bk_tree_search(BKTree *tree, const char *query, int threshold,
                   SearchResult *results, int max_results)
{
    int count = 0;
    tree->nodes_visited = 0;

    bk_search_recursive(
        tree->root, query, (int)strlen(query),
        threshold, results, &count, max_results,
        &tree->nodes_visited
    );

    return count;
}

/* ---------- 统计信息 ---------- */

void bk_tree_print_stats(BKTree *tree)
{
    printf("BK-tree 统计:\n");
    printf("  词条数量:       %d\n", tree->size);
    printf("  上次查询访问:   %d 个节点 (%.1f%%)\n",
           tree->nodes_visited,
           tree->size > 0 ? 100.0 * tree->nodes_visited / tree->size : 0);
}

/* ---------- 演示 ---------- */

int demo_bk_tree(void)
{
    /* 构建一个小词典 */
    const char *dict[] = {
        "apple", "apply", "ape", "able", "april",
        "orange", "range", "organ", "origin",
        "python", "typhon", "pylon", "prion",
        "search", "starch", "scorch", "parch",
        "hello", "hallo", "hullo", "help",
        "world", "word", "worm", "work", "worst",
        "fuzzy", "fussy", "furry", "funny", "fully",
        "spell", "shell", "smell", "swell", "skill",
        "check", "chess", "chest", "cheek", "creek",
        NULL
    };

    BKTree *tree = bk_tree_create();
    for (int i = 0; dict[i]; i++)
        bk_tree_insert(tree, dict[i]);

    printf("词典加载完成，共 %d 个词\n\n", tree->size);

    /* 模糊搜索测试 */
    const char *queries[] = {"pytohn", "serch", "wrold", "helo", "aple"};
    int num_queries = 5;

    for (int q = 0; q < num_queries; q++) {
        SearchResult results[MAX_RESULTS];
        int count = bk_tree_search(tree, queries[q], 2, results, MAX_RESULTS);

        printf("查询 \"%s\" (阈值=2):\n", queries[q]);
        for (int i = 0; i < count; i++)
            printf("  距离=%d: %s\n", results[i].distance, results[i].word);
        bk_tree_print_stats(tree);
        printf("\n");
    }

    bk_tree_destroy(tree);
    return 0;
}

BK-tree 的局限性

BK-tree 并非完美。它有几个需要注意的问题：

1. 构建不平衡：树的形状高度依赖于插入顺序和根节点的选择。随机选根通常比选第一个词更好。
2. 缓存不友好：树的节点分散在堆上，指针追踪（pointer chasing）导致大量缓存缺失。
3. 高阈值退化：当 k >= 3 时，剪枝效果急剧下降，可能退化到接近暴力搜索。
4. 不适合并行：树结构的递归搜索难以有效并行化。

对于大规模词典（百万级以上），更好的选择可能是 VP-tree、倒排索引、或下面要讲的 Levenshtein 自动机。

七、Levenshtein 自动机：NFA 到 DFA 的快速匹配

Levenshtein 自动机是另一种截然不同的模糊匹配方法。它的核心思想是：给定模式串 p 和距离阈值 k，构造一个自动机 A(p, k)，使得 A 能够接受且仅接受与 p 的编辑距离 <= k 的所有字符串。

NFA 构造

Levenshtein NFA 的状态空间是一个 (m+1) x (k+1) 的网格，状态 (i, e) 表示”已经处理了模式串的前 i 个字符，消耗了 e 次编辑操作”。

转移规则： 1. 匹配：如果当前输入字符 c == p[i]，从 (i, e) 转移到 (i+1, e)。 2. 替换：从 (i, e) 通过任意字符转移到 (i+1, e+1)。 3. 插入：从 (i, e) 通过任意字符转移到 (i, e+1)（模式串位置不变）。 4. 删除：从 (i, e) 通过 epsilon（不消耗输入）转移到 (i+1, e+1)。

初始状态是 (0, 0)。接受状态是所有 (m, e) 其中 e <= k，以及通过”尾部通配”到达的状态。

从 NFA 到 DFA

NFA 直接使用效率不高（需要追踪状态集合），所以通常会通过子集构造法（subset construction）将其转换为 DFA。转换后的 DFA 状态数是有限的，而且 Schulz 和 Mihov 在 2002 年证明了一个重要结论：对于固定的 k，Levenshtein DFA 的状态模式是可以参数化的（parametric states），不依赖于具体的模式串内容，只依赖于模式串长度。

这意味着可以预计算通用的转移表，然后对任何具体的模式串在 O(m) 时间内实例化 DFA。

与 Trie 结合

Levenshtein 自动机最强大的用法是与字典的 Trie 结构结合。将自动机和 Trie 同时遍历，自动机的状态跟随 Trie 的分支前进。当自动机进入死状态（距离必然超过 k）时，整棵子树都可以剪掉。

这种方法的查询时间只依赖于自动机的大小和结果集的大小，与词典的规模几乎无关。对于 k=2 的查询，在百万级词典上通常可以在微秒级完成。

伪代码

function levenshtein_automaton_search(trie_root, pattern, k):
    nfa = build_levenshtein_nfa(pattern, k)
    dfa = subset_construction(nfa)

    results = []
    stack = [(trie_root, dfa.start_state, "")]

    while stack is not empty:
        (trie_node, dfa_state, prefix) = stack.pop()

        if dfa_state is accepting:
            results.append(prefix + remaining_chars(trie_node))

        for each child (char c, child_node) of trie_node:
            next_state = dfa.transition(dfa_state, c)
            if next_state is not dead:
                stack.push((child_node, next_state, prefix + c))

    return results

实际应用

Lucene（Elasticsearch 的底层库）从 4.0 版本开始就使用 Levenshtein 自动机来实现 fuzzy query。它预编译了 k=1 和 k=2 的参数化 DFA 转移表，然后在查询时与倒排索引的 term dictionary（基于 FST）进行交叉遍历。这使得 Elasticsearch 的模糊查询可以在毫秒级返回结果，即使面对数百万个不同的 term。

八、基准测试与性能对比

理论分析需要实验验证。以下是在一台普通开发机器（Intel i7-12700K，DDR5-4800，GCC 12 -O2）上的基准测试结果。

测试一：单次距离计算

测试两个字符串之间编辑距离的计算速度。模式串和文本串从英文词典中随机选取。

Myers 算法在 m<=64 时展现出压倒性的优势。对于 m=64 的情况，它比朴素 DP 快了 100 倍以上。

测试二：词典模糊搜索

词典大小 50 万个英文单词，查询 1000 个随机拼写错误的词，阈值 k=2。

几个关键观察：

1. 距离函数的选择至关重要。在暴力搜索中，仅将距离函数从朴素 DP 换成 Myers，就获得了 10 倍加速。
2. BK-tree 的剪枝效果显著。它将比较次数从 50 万降到了 7.5 万（约 15%），再配合 Myers 算法，查询时间降到 1.8 毫秒。
3. Levenshtein 自动机 + Trie 是量级碾压。0.05 毫秒的查询时间意味着每秒可以处理两万次模糊查询，完全满足在线服务的需求。

测试三：不同阈值 k 的影响

使用 BK-tree + Myers，词典大小 50 万：

可以清楚地看到，k 每增加 1，搜索范围大致扩大 3-5 倍。这就是为什么实际系统中 k 很少超过 2 的原因。

九、实战应用场景

编辑距离和模糊匹配在工业界有广泛的应用。以下是几个最重要的场景。

Elasticsearch fuzzy query

Elasticsearch 通过 fuzziness 参数支持模糊查询：

{
  "query": {
    "match": {
      "title": {
        "query": "pytohn",
        "fuzziness": "AUTO"
      }
    }
  }
}

fuzziness: "AUTO" 的规则是： - 字符串长度 0-2：必须精确匹配（k=0） - 字符串长度 3-5：允许 1 次编辑（k=1） - 字符串长度 >5：允许 2 次编辑（k=2）

底层实现使用 Levenshtein 自动机与 FST（Finite State Transducer）结构的 term dictionary 进行交叉遍历。注意 Elasticsearch 使用的是 Damerau-Levenshtein 距离（包含转置），这对打字错误更友好。

拼写检查器

现代拼写检查器的典型架构：

用户输入 "recieve"
    |
    v
词典查找 --> 未命中
    |
    v
生成候选词（编辑距离 <= 2）
    |
    v
候选排序（语言模型概率 + 距离权重）
    |
    v
返回建议："receive", "relieve", "retrieve"

Peter Norvig 的经典文章 “How to Write a Spelling Corrector” 用 Python 21 行代码实现了一个惊人有效的拼写纠错器。其核心思想是：

1. 生成所有编辑距离为 1 的候选词（对于长度为 n 的词，约有 54n+25 个候选）。
2. 在词典中过滤存在的候选词。
3. 如果没有找到，继续生成距离为 2 的候选词。
4. 按词频排序返回最可能的纠正。

这种”生成-过滤”的方法对小阈值非常有效，但随着 k 的增大，候选词数量指数增长，此时 BK-tree 或自动机方法更有优势。

diff 算法

Unix diff 和 Git diff 的核心是寻找两个文件的最长公共子序列（LCS）或等价地寻找最小编辑距离。Myers 在 1986 年发表的 O(ND) 算法是这个领域的里程碑：

/* Myers 差异算法的核心思想 */
/*
 * 在编辑图上寻找从 (0,0) 到 (N,M) 的最短路径
 * 对角线移动（匹配）是免费的
 * 水平移动（删除）和垂直移动（插入）代价为 1
 * D 是编辑距离，搜索按 D=0, 1, 2, ... 逐层展开
 */

Git 默认使用这个算法，当文件差异较小时效率极高（时间复杂度 O((N+M)D)，D 是差异大小）。

DNA 序列比对

在生物信息学中，编辑距离的变体被广泛用于 DNA 和蛋白质序列比对。与字符串拼写纠错不同，生物序列比对通常使用加权的编辑距离（不同的替换有不同的代价，通过替换矩阵如 BLOSUM62 定义），并且需要处理更长的序列（数千到数百万碱基）。

BLAST、BWA、Bowtie 等比对工具都在内部使用了编辑距离或其变体。其中 Myers 的位并行算法在短读段（short reads）比对中尤为常用。

模糊命令行补全

很多现代命令行工具使用编辑距离来实现”你是不是想输入…“的功能：

$ git comit
git: 'comit' is not a git command. Did you mean 'commit'?

$ cargo buidl
error: no such command: `buidl`

    Did you mean `build`?

这些通常使用简单的有界编辑距离实现，阈值设为 1 或 2，在几十个候选命令中搜索即可。

十、工程实践中的陷阱与技巧

在将编辑距离算法部署到生产环境时，有许多容易被忽视的细节。以下是我在实践中总结的经验。

常见陷阱一览

性能优化清单

1. 选择正确的算法：短字符串（<= 64）用 Myers，长字符串用空间优化 DP。
2. 提前终止：如果只关心距离是否 <= k，使用有界版本。
3. 长度差剪枝：如果 |m - n| > k，距离必然 > k，直接跳过。
4. 公共前缀/后缀剥离：两个串的公共前缀和后缀不影响编辑距离，可以先去掉再计算。
5. SIMD 向量化：对于批量计算，可以用 SSE/AVX 同时处理多个距离计算。
6. 缓存预热：对于 BK-tree，将高频查询对应的子树预加载到 CPU 缓存中。

公共前缀后缀剥离

这是一个被严重低估的优化。实现简单，效果显著：

int levenshtein_stripped(const char *s, int m, const char *t, int n)
{
    /* 剥离公共前缀 */
    while (m > 0 && n > 0 && *s == *t) {
        s++; t++; m--; n--;
    }

    /* 剥离公共后缀 */
    while (m > 0 && n > 0 && s[m - 1] == t[n - 1]) {
        m--; n--;
    }

    /* 对剥离后的子串计算距离 */
    if (m == 0) return n;
    if (n == 0) return m;

    return levenshtein_one_row(s, m, t, n);
}

对于 “internationalization” 和 “internationalisation”（s/z 之差），公共前缀 “internationali” 和公共后缀 “ation” 被剥离后，只需要计算 “z” 和 “s” 的距离（= 1）。原本需要 20x21 = 420 个 DP 单元格，优化后只需要 1x2 = 2 个。

十一、延伸话题与进阶阅读

加权编辑距离

标准 Levenshtein 距离假设所有操作的代价相同，但在实际应用中这往往不够准确。例如：

• 键盘上相邻按键的替换（如 a->s）应该比远距离替换（如 a->p）代价更低。
• 元音之间的替换（如 e->i）在某些语言中比辅音替换更常见。
• 删除重复字符（如 “helllo”->“hello”）的代价应该更低。

加权编辑距离只需要修改 DP 的转移方程，将固定的代价 1 替换为查表得到的代价值。Wagner-Fischer 算法完全适用，只是状态转移中的 +1 变成了 +w(op)。

编辑距离的近似算法

对于超长序列（如完整的基因组比对），即使 O(mn) 的精确算法也太慢了。研究者提出了一些近似和随机化方法：

• 带状 DP（Banded DP）：只计算主对角线附近宽度为 2k+1 的带状区域，时间复杂度降为 O(nk)。
• 嵌入方法：将字符串嵌入到低维空间，用欧几里得距离近似编辑距离。
• 采样方法：随机采样子串进行比较，估算全局距离。

这些方法在精度和速度之间做出权衡，适用于不同的场景。

树编辑距离

编辑距离的概念可以从字符串推广到树结构。树编辑距离度量将一棵有标签的树变换为另一棵树所需的最少节点操作次数（插入、删除、重标记）。这在 XML/HTML 比较、语法树分析、RNA 二级结构比对等领域有重要应用。

Zhang-Shasha 算法（1989）可以在 O(n1 * n2 * min(d1, l1) * min(d2, l2)) 时间内计算树编辑距离，其中 n 是节点数，d 是深度，l 是叶子数。

推荐论文与资源

1. Levenshtein (1965). “Binary codes capable of correcting deletions, insertions, and reversals.” 开山之作。
2. Wagner & Fischer (1974). “The string-to-string correction problem.” DP 算法的奠基论文。
3. Myers (1999). “A fast bit-vector algorithm for approximate string matching.” 位并行算法。
4. Myers (1986). “An O(ND) difference algorithm and its variations.” Git diff 的基础。
5. Schulz & Mihov (2002). “Fast string correction with Levenshtein automata.” 参数化自动机。
6. Norvig (2007). “How to Write a Spelling Corrector.” 经典教程。
7. Burkhard & Keller (1973). “Some approaches to best-match file searching.” BK-tree 原论文。

十二、个人思考

写到这里，回顾一下编辑距离这个话题，有几点感触。

关于算法选择的务实态度。在大多数拼写纠错的实际场景中，你要处理的模式串不超过 20 个字符，阈值不超过 2。这种情况下，一个经过仔细优化的有界 DP（加上公共前缀后缀剥离和长度差剪枝）已经足够快了。Myers 算法虽然理论上更优，但实现复杂度也更高。不要为了追求理论最优而忽视工程上的简单性。当然，如果你是在写一个通用库、一个搜索引擎的底层组件，那么投入时间实现 Myers 是值得的。

关于度量空间的深层意义。Levenshtein 距离满足度量公理这件事情看似平凡，实则深刻。它意味着整个度量空间理论的工具箱（BK-tree、VP-tree、cover tree、各种近邻搜索算法）都可以直接拿来用。如果你定义了一个不满足三角不等式的”距离函数”，你就失去了这一切。这也是为什么 OSA 距离（不满足三角不等式）在实际中使用受限的原因。在设计距离函数时，三角不等式是一条值得守护的底线。

关于自动机方法的优雅。Levenshtein 自动机把一个看似需要动态计算的问题转化为静态的状态机匹配。这种”编译”思维在计算机科学中屡见不鲜：正则表达式编译为 DFA，SQL 查询编译为执行计划，着色器编译为 GPU 指令。预计算和编译永远是性能优化的最强武器之一。

关于中文编辑距离的特殊性。中文的编辑距离比英文复杂得多。首先是字符粒度的问题：按字符还是按拼音？“北京” 和 “被禁” 字面上完全不同（编辑距离 = 2），但拼音几乎相同（“beijing” vs “beijing”，距离 = 0）。其次是形近字的问题：“己” “已” “巳” 字形相似但 Unicode 码点完全不同。实际的中文拼写纠错系统通常需要综合考虑拼音距离、字形距离和语义距离，远比纯英文系统复杂。

关于”足够好”的工程哲学。拼写纠错不需要完美。用户搜索 “pytohn”，你给出 “python” 作为建议就够了，不需要穷举所有编辑距离 <= 2 的候选词。在工程实践中，“前 5 个最可能的纠正”比”所有可能的纠正”重要得多。这意味着你可以在搜索过程中激进地剪枝，在候选排序中引入语言模型和上下文信息。完美的编辑距离计算只是整个纠错管道中的一个环节，而且往往不是瓶颈所在。

最后，如果你正在开发一个需要模糊匹配的系统，我的建议是这样的优先级：

1. 先评估数据规模和延迟要求。
2. 如果词典 < 10 万且 k <= 2，用暴力 + Myers 可能就够了。
3. 如果词典更大，考虑 BK-tree 或 Trie + 自动机。
4. 如果在已有的搜索引擎上，直接用 Elasticsearch 的 fuzzy query。
5. 无论如何，记得做大小写归一化和 Unicode 规范化。

这些算法和数据结构发明于上世纪六七十年代，但至今仍在驱动着我们每天使用的搜索引擎和文本编辑器。它们是那种”不会过时的知识”，值得花时间深入理解。

上一篇: BWT 与 FM-index

下一篇: 字符串哈希

相关阅读: - AC 自动机 - 区间 DP