汚れたQRコードがなぜ読めるのか、傷ついたCDがなぜ再生できるのか、そして地球から数億キロ離れた探査機の画像がなぜ破綻せずに届くのか。答えはどれも同じで、リードソロモン符号（RS符号）という誤り訂正符号が背後で働いているからです。ところが「RS符号は $t$ 個までのシンボル誤りを訂正できる」という説明を聞いたとき、多くの人が同じところで引っかかります。訂正できることは分かった。でも、255個のシンボルのうちどの8個が壊れているかを、たった16個の数値（シンドローム）からどうやって特定するのか？
素朴に考えると、255個から8個を選ぶ組み合わせは約 $4 \times 10^{14}$ 通りあります。総当たりは論外です。しかも「どの位置が壊れたか」だけでなく「どれだけ壊れたか（誤り値）」も同時に求めなければなりません。未知数は $8 + 8 = 16$ 個、手がかりも16個。数の上では釣り合っていますが、方程式が未知数の冪乗を含む非線形連立方程式になるため、そのままでは解けません。
この難問を、「与えられた数列を生成する最短の線形フィードバックシフトレジスタ（LFSR）を求めよ」という、まったく別の問題にすり替えて一気に解いてしまうのが Berlekamp-Masseyアルゴリズム（BMアルゴリズム）です。1968年にBerlekampがBCH符号の復号法として発表し、翌1969年にMasseyが「これはLFSRの合成問題そのものだ」と読み替えて、驚くほど短いアルゴリズムに整理しました。計算量は $O(t^2)$ で、$t=8$ なら数十回の有限体乗算で終わります。
BMアルゴリズムを理解すると、次のような場面がつながって見えてきます。
- ストレージと放送: CD・DVD・Blu-ray、QRコード、DVB放送、SSDのECC。これらはすべてRS符号を使っており、その復号器の心臓部がBMアルゴリズムです
- 深宇宙通信: ボイジャーやカッシーニをはじめとするCCSDS規格の宇宙機通信では、RS(255,223) が畳み込み符号と連結して使われてきました。低SNR環境でバースト誤りを叩き落とす役割を担っています
- 暗号解析: LFSRを鍵ストリーム生成器に使うストリーム暗号の安全性尺度である「線形複雑度」は、まさにBMアルゴリズムが返す $L$ です。BMで短いLFSRが見つかる鍵ストリームは、そこから先が全部予測できてしまいます
- 有限体上の高速アルゴリズム: Padé近似・部分分数展開・拡張ユークリッド互除法と数学的に等価な構造を持ち、記号計算の分野でも使われます
本記事の内容
- 誤りパターンとシンドロームの関係 — なぜ素朴な連立方程式では解けないのか
- 誤り位置多項式 $\Lambda(x)$ と鍵方程式 $\Lambda(x)S(x) \equiv \Omega(x) \pmod{x^{2t}}$ の導出（省略なし）
- 鍵方程式が「シンドローム列を生成する最短LFSRを求める問題」と等価であることの証明
- ディスクレパンシー $d$、接続多項式の更新 $\Lambda^{(n)} = \Lambda^{(n-1)} – (d/b)\,x^m B(x)$、長さ更新則 $L \leftarrow n+1-L$ の意味と最短性の帰納法による証明
- GF(16) 上の RS(15,11) による完全な手計算例
- GF($2^8$) の四則演算・RS(255,239)の符号化・BM・Chien探索・Forney公式のPythonフルスクラッチ実装
- 誤り数を0〜12個まで振ったときに、訂正成功率が $t=8$ を境に崖のように崩れることの実測
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
- リードソロモン符号の理論 — バースト誤りに強い符号を数学から理解する — 本記事が復号する対象そのもの。生成多項式とガロア体上の符号語の定義
- BCH符号 — ガロア体上で設計する強力な誤り訂正符号 — BM復号はBCH符号でも同じ形で使えます。BCH限界と設計距離の考え方
- ハミング符号とは？誤り訂正の仕組みを図解でわかりやすく解説 — シンドロームという概念の一番やさしい入口
- 巡回符号の理論 — 生成多項式とシフトレジスタ符号化を基礎から解説 — 多項式で符号を扱う流儀と、シフトレジスタによる実装
復号の全体像 — 4つのステージ
具体的な数式に入る前に、RS符号の復号がどんな流れになっているかを俯瞰しておきましょう。復号器は次の4つのステージを順に通ります。
この図が本記事の地図です。受信語から訂正済み符号語までの流れが4つの箱に分かれていて、そのうち赤で囲った2番目だけが「代入するだけ」では済まない難所になっています。シンドローム計算もChien探索もForney公式も、やっていることは有限体の元を多項式に代入する単純作業で、実装も理論も難しくありません。本記事の残りの大半は、この赤い箱の中身、すなわち「16個の数値から $\Lambda(x)$ という多項式を復元する」ためだけに費やされます。
- シンドローム計算: 受信語 $r(x)$ を生成多項式の根 $\alpha, \alpha^2, \dots, \alpha^{2t}$ に代入して $2t$ 個の値 $S_1, \dots, S_{2t}$ を得る。無誤りなら全部ゼロ
- 誤り位置多項式の決定: シンドローム列から $\Lambda(x)$ を求める。ここがBMアルゴリズムの担当
- Chien探索: $\Lambda(x)$ の根を全数探索で見つけ、誤りが起きた位置に翻訳する
- Forney公式: 誤り値を閉じた式で計算し、受信語から差し引いて訂正完了
面白いのは、この4ステージのうち3つは「代入して評価するだけ」の単純作業だという点です。シンドロームは代入、Chien探索も代入、Forneyも代入。難しいのはステージ2、つまり16個の数値から $\Lambda(x)$ という多項式を復元する部分に集中しています。BMアルゴリズムはこの一点だけを担当する、いわば復号器の頭脳です。
例えるなら、犯人捜しに似ています。シンドロームは「現場に残された16個の痕跡」です。Chien探索は「容疑者リストを1人ずつ当たる」地道な作業、Forneyは「犯人が確定したあとで被害額を計算する」経理作業。BMは、16個の痕跡から犯人の人数と、犯人グループを特徴づける方程式を一気に割り出す推理パートに当たります。
まずは痕跡がどうやって残るのか、つまりシンドロームと誤りの関係から見ていきましょう。
誤りの数学 — シンドロームは何を記録しているのか
誤りパターンを多項式で書く
RS符号では、符号長 $n$ の符号語をGF($q$)上の多項式として扱います。以降、$q = 2^m$ とし、$\alpha$ をGF($2^m$)の原始元とします。符号語多項式 $c(x)$ は生成多項式
$$ g(x) = \prod_{j=1}^{2t} (x – \alpha^j) $$
の倍数として定義されます。つまり符号語は $c(\alpha) = c(\alpha^2) = \cdots = c(\alpha^{2t}) = 0$ を満たす多項式です。これが「符号語であることの証明書」になっています。
通信路を通ると、受信語は
$$ r(x) = c(x) + e(x) $$
となります。ここで $e(x)$ が誤りパターンです。いま $\nu$ 個のシンボルが壊れたとして、壊れた位置の指数を $i_1, i_2, \dots, i_\nu$、そこに乗った誤りの大きさを $e_{i_1}, \dots, e_{i_\nu}$（いずれも非零）と書くと
$$ e(x) = \sum_{k=1}^{\nu} e_{i_k} x^{i_k} $$
です。ここで記号を整理します。誤り位置数（error locator）を $X_k = \alpha^{i_k}$、誤り値（error value）を $Y_k = e_{i_k}$ と置きます。位置という「整数」を、$\alpha$ の冪という「有限体の元」に変換したのがポイントです。整数の $i_k$ は扱いにくいですが、$X_k = \alpha^{i_k}$ なら体の中で掛け算や割り算ができます。
左の(a)は、後の節で手計算するGF(16)の例です。位置 $i_1=10$ と $i_2=3$ に誤りが加わっています（$i=3$ は元のパリティ値 $\alpha^5$ と同じ値が加わって偶然ゼロになっていますが、「壊れている」ことに変わりはありません）。右の(b)は、その位置を $X_k = \alpha^{i_k}$ という体の元に翻訳したものです。GF(16)の非零元は $\alpha^0$ から $\alpha^{14}$ までの15個が円環をなしていて、位置番号 $i$ と体の元 $\alpha^i$ が一対一に対応します。この翻訳のおかげで、以降は「位置」を掛けたり割ったりできるようになり、代数的な操作の対象になります。
シンドロームは冪和である
シンドロームは受信語を根に代入した値です。
$$ S_j = r(\alpha^j) = c(\alpha^j) + e(\alpha^j) = 0 + e(\alpha^j) = \sum_{k=1}^{\nu} e_{i_k} (\alpha^j)^{i_k} $$
$c(\alpha^j) = 0$ が効いて、シンドロームには誤りの情報だけが残ります。
上段の(a)(b)を見ると、255シンボルもの中身がぎっしり詰まった符号語であっても、そのシンドローム16個は完全にゼロです。ユーザーデータがどれだけ複雑でも、$g(x)$ の根で評価すれば必ず消えるように符号化してあるからです。下段の(c)(d)では2シンボルだけを壊しました。受信語の見た目はほとんど変わりませんが、シンドロームは一気に非ゼロの値で埋まります。255個の情報のうち、誤りに関する成分だけが16個の数値に圧縮されて抽出される — これが誤り訂正符号の設計思想の核心です。
つまり符号語の中身（ユーザーデータ）は完全に消え去り、誤りの痕跡だけが抽出されるのです。$(\alpha^j)^{i_k} = (\alpha^{i_k})^j = X_k^j$ と書き換えると、非常にきれいな形になります。
$$ \begin{equation} S_j = \sum_{k=1}^{\nu} Y_k X_k^{\,j}, \qquad j = 1, 2, \dots, 2t \end{equation} $$
これは統計でいう「冪和」（power sum）の形です。$X_k$ を「$k$ 番目の犯人のID」、$Y_k$ を「その犯人の寄与の重み」と見れば、$S_1$ は重み付き1乗和、$S_2$ は重み付き2乗和……というふうに、同じ犯人グループを次々と違う角度から見た統計量が $2t$ 個並んでいることになります。
なぜこのままでは解けないのか
式(1)を未知数 $X_k, Y_k$ に関する連立方程式と見ると、$\nu \le t$ のとき未知数は $2\nu \le 2t$ 個、方程式は $2t$ 個です。数の上では足りています。ところが問題があります。$Y_k$ については確かに線形ですが、$X_k$ については $X_k^j$ という高次の冪が現れる非線形方程式なのです。ガウス消去法のような線形代数の道具がそのままでは使えません。
しかも、$\nu$（誤りの個数）そのものが未知です。方程式の本数は決まっていても、未知数の個数が未知という状況です。誤りが3個なのか5個なのか分からないまま、$\{X_k\}$ と $\{Y_k\}$ を同時に決めなければなりません。
この二重の困難を突破するために、先人は発想を転換しました。$X_k$ を直接求めるのをあきらめ、$X_k$ たちを根に持つ多項式の係数を求めることにしたのです。係数なら線形方程式で決まるかもしれない、という賭けです。次節でこの賭けが見事に当たることを見ていきます。
誤り位置多項式と鍵方程式
誤り位置多項式の定義
$n$ 次方程式の解を直接求めるのは難しくても、解を根に持つ多項式の係数（基本対称式）なら扱いやすい。これは高校で習う解と係数の関係と同じ発想です。$x^2 – 5x + 6 = 0$ の解が2と3であることを求めるより、「和が5、積が6」という情報を先に手に入れるほうが簡単な場面があります。
そこで、誤り位置数の逆数を根に持つ多項式を定義します。
$$ \begin{equation} \Lambda(x) = \prod_{k=1}^{\nu} (1 – X_k x) = 1 + \Lambda_1 x + \Lambda_2 x^2 + \cdots + \Lambda_\nu x^\nu \end{equation} $$
これを誤り位置多項式（error locator polynomial）と呼びます。定義から次の3つの性質が読み取れます。
- 次数が誤り個数を教える: $\deg \Lambda = \nu$ です。$\Lambda(x)$ が求まれば、誤りが何個あったかが自動的に分かります
- 定数項が1: $\Lambda(0) = 1$ です。これは正規化であり、後で「更新の起点」として効いてきます
- 根が誤り位置を教える: $\Lambda(X_k^{-1}) = 0$ です。逆数を根にしたのは、定数項を1に固定するためです（$\prod(x – X_k)$ だと定数項が $\prod X_k$ になってしまう）
左の(a)は、2つの誤り位置数 $X_1=\alpha^{10}$、$X_2=\alpha^3$ から $\Lambda(x)$ を組み立てる流れです。展開した結果 $1+\alpha^{12}x+\alpha^{13}x^2$ には、誤り位置の情報が「係数」という形で埋め込まれています。右の(b)は、その $\Lambda(x)$ に候補 $\alpha^{-i}$（$i=0,\dots,14$）を総当たりで代入した値です。15個のうち $i=3$ と $i=10$ のちょうど2つだけがゼロになり、実際に壊れた位置とぴたりと一致しています。係数さえ手に入れば、あとはこの棒グラフでゼロを探すだけで誤り位置が確定するわけです。
つまり $\Lambda(x)$ さえ手に入れば、根を探すだけで誤り位置がすべて分かります。有限体の元は $2^m – 1$ 個しかないので、全数代入（＝Chien探索）が現実的な計算量で回ります。問題は「$\Lambda(x)$ をシンドロームからどう作るか」です。
シンドローム多項式の閉じた形
シンドローム $2t$ 個を係数に持つ多項式を定義します。
$$ S(x) = S_1 + S_2 x + S_3 x^2 + \cdots + S_{2t} x^{2t-1} = \sum_{j=1}^{2t} S_j x^{j-1} $$
ここに式(1)を代入して、和の順序を入れ替えます。
$$ S(x) = \sum_{j=1}^{2t} \left( \sum_{k=1}^{\nu} Y_k X_k^{\,j} \right) x^{j-1} = \sum_{k=1}^{\nu} Y_k \sum_{j=1}^{2t} X_k^{\,j} x^{j-1} $$
内側の和で $X_k$ を1つくくり出すと、$\sum_{j=1}^{2t} X_k^{\,j} x^{j-1} = X_k \sum_{j=1}^{2t} (X_k x)^{j-1}$ となります。これは公比 $X_k x$、初項1、項数 $2t$ の等比数列の和ですから、等比級数の公式が使えます。
$$ \sum_{j=1}^{2t} (X_k x)^{j-1} = \frac{1 – (X_k x)^{2t}}{1 – X_k x} $$
これを戻すと、シンドローム多項式の閉じた形が得られます。
$$ \begin{equation} S(x) = \sum_{k=1}^{\nu} \frac{Y_k X_k \left( 1 – (X_k x)^{2t} \right)}{1 – X_k x} \end{equation} $$
分母に $1 – X_k x$ が現れました。これは $\Lambda(x)$ の因子そのものです。ここで「$\Lambda(x)$ を掛ければ分母が払える」という道筋が見えます。
鍵方程式の導出
式(3)の両辺に $\Lambda(x) = \prod_{l=1}^{\nu}(1 – X_l x)$ を掛けます。$k$ 番目の項では分母 $1 – X_k x$ が約分され、$k$ 以外の因子だけが残ります。
$$ \Lambda(x) S(x) = \sum_{k=1}^{\nu} Y_k X_k \left( 1 – (X_k x)^{2t} \right) \prod_{l \ne k} (1 – X_l x) $$
ここで $(X_k x)^{2t} = X_k^{2t} x^{2t}$ という項に注目します。この項は必ず $x^{2t}$ の因子を持ちますから、$x^{2t}$ を法として考えれば消えてしまいます。つまり $\bmod\ x^{2t}$ を取ると
$$ \Lambda(x) S(x) \equiv \sum_{k=1}^{\nu} Y_k X_k \prod_{l \ne k} (1 – X_l x) \pmod{x^{2t}} $$
となります。右辺を誤り評価多項式（error evaluator polynomial）と名づけ、$\Omega(x)$ と書きます。
$$ \begin{equation} \Omega(x) = \sum_{k=1}^{\nu} Y_k X_k \prod_{l \ne k} (1 – X_l x) \end{equation} $$
以上で、RS/BCH復号の中心にある関係式が得られました。
$$ \begin{equation} \boxed{\ \Lambda(x) S(x) \equiv \Omega(x) \pmod{x^{2t}}\ } \end{equation} $$
これを鍵方程式（key equation）と呼びます。既知はシンドローム多項式 $S(x)$ のみ、未知は $\Lambda(x)$ と $\Omega(x)$ の2つです。
鍵方程式が「解ける」理由
未知が2つあるのに解けるのは、次数に強い制約があるからです。式(2)より $\deg \Lambda = \nu$。式(4)を見ると、各項は $\nu – 1$ 個の1次因子の積なので $\deg \Omega \le \nu – 1$ です。まとめると
$$ \deg \Omega \le \nu – 1 < \nu = \deg \Lambda \le t $$
という不等式が成り立ちます。「$\Omega$ の次数が $\Lambda$ より真に小さい」「どちらも $t$ 以下」という2条件を満たす解は、実は本質的に一意であることが証明できます（$\Lambda$ と $\Omega$ が互いに素であるという条件を加えれば完全に一意）。だから鍵方程式は well-posed な問題なのです。
もうひとつ重要な性質があります。式(4)で $x = X_k^{-1}$ を代入すると、$l \ne k$ の因子はどれもゼロにならず、$l = k$ の項だけが生き残ります。
$$ \Omega(X_k^{-1}) = Y_k X_k \prod_{l \ne k} (1 – X_l X_k^{-1}) \ne 0 $$
つまり $\Omega$ は $\Lambda$ の根で消えない。この事実が後のForney公式（誤り値の閉形式）を支えます。
ここまでで「$\Lambda(x)$ を求めれば勝ち」という構図と、その $\Lambda(x)$ が満たす方程式が分かりました。しかし式(5)はまだ「合同式」であって、アルゴリズムではありません。次節で、これを逐次的に解けるかたちに翻訳します。
鍵方程式はLFSRの合成問題である
合同式を係数ごとに書き下す
鍵方程式(5)の意味は「左辺と右辺の $x^0$ から $x^{2t-1}$ までの係数が全部一致する」ことです。そこで係数を比較してみましょう。左辺 $\Lambda(x)S(x)$ の $x^{n}$ の係数は、畳み込みですから
$$ [x^n]\,\Lambda(x)S(x) = \sum_{i=0}^{\min(n,\nu)} \Lambda_i \, S_{n+1-i} $$
です（$\Lambda_0 = 1$、$S_j$ の添字は1始まりなので $x^{n-i}$ の係数は $S_{n-i+1}$）。一方、右辺の $\Omega(x)$ は $\deg \Omega \le \nu – 1$ なので、$n \ge \nu$ の範囲では右辺の係数はすべてゼロです。したがって $n = \nu, \nu+1, \dots, 2t-1$ に対して
$$ \sum_{i=0}^{\nu} \Lambda_i \, S_{n+1-i} = 0 $$
が成り立ちます。
この帯図は鍵方程式を「係数の一覧表」として並べ直したものです（見やすさのため $2t=8$、$\nu=2$ の場合）。上段が $\Lambda(x)S(x)$ の係数、中段が $\Omega(x)$ の係数、下段がその比較結果です。ポイントは2本の壁です。紫の壁（$x^{2t}$）より右は合同式の法によって切り捨てられるので、そもそも情報がありません。赤の壁（$\deg\Omega \le \nu-1$）より左は $\Omega$ の係数と等しくなりますが、$\Omega$ 自体が未知なので何も言えません。使えるのは2本の壁に挟まれた中央の領域だけで、そこでは右辺がゼロだと確定している — この $2t-\nu$ 本の「ゼロである」という式こそが、次に見る漸化式の正体です。
$i = 0$ の項（$\Lambda_0 = 1$）を分離し、添字を $n+1 \to n$ とずらして書き直すと、次の漸化式になります。
$$ \begin{equation} S_n = -\sum_{i=1}^{\nu} \Lambda_i \, S_{n-i}, \qquad n = \nu+1, \nu+2, \dots, 2t \end{equation} $$
標数2の体（GF($2^m$)）では $-1 = +1$ ですから、実装上は符号を気にせず
$$ S_n = \sum_{i=1}^{\nu} \Lambda_i S_{n-i} \quad (\text{GF}(2^m) \text{ の場合}) $$
と書けます。
これはLFSRそのもの
式(6)をじっと見てください。「数列の次の項が、直前 $\nu$ 項の一次結合で決まる」という形です。これは線形フィードバックシフトレジスタ（LFSR）の動作そのものです。
LFSRとは、$\nu$ 段のレジスタに値を蓄え、タップ係数 $\Lambda_1, \dots, \Lambda_\nu$ で重み付けした和をフィードバックして、次のシンボルを吐き出す回路です。$\Lambda(x)$ はこの回路の配線図であり、この文脈では接続多項式（connection polynomial）と呼ばれます。$L = \nu$ がレジスタの段数（LFSRの長さ）です。
回路図として描くと、$\Lambda(x)$ の係数が何であるかが一目で分かります。青い箱が直前4個のシンボルを覚えているレジスタ、オレンジの丸が係数 $\Lambda_1,\dots,\Lambda_4$ を掛ける乗算器、紫の丸がそれらを足し合わせる加算器（GF($2^m$) なのでXOR）です。$\Lambda(x)$ の係数はそのまま「どのタップにいくつを掛けるか」という配線の設定値になっています。したがって「$\Lambda(x)$ を求めよ」という問題は「この回路の配線をどう設定すれば、与えられたシンドローム列を吐き出すか」という設計問題と同じものです。しかも段数 $L$ は誤り個数 $\nu$ と一致するので、最短の回路を見つけることが誤り個数を当てることにもなります。
つまり、鍵方程式を解くことは次の問題と等価です。
問題: 数列 $S_1, S_2, \dots, S_{2t}$ が与えられたとき、この数列を生成できるLFSRのうち、最短のもの（＝ $L$ が最小のもの）の接続多項式 $\Lambda(x)$ と長さ $L$ を求めよ。
これがMasseyの読み替えです。「誤り位置を探す」という符号理論の問題が、「数列を再現する最小の回路を設計する」というシステム同定の問題になりました。
なぜ「最短」でなければならないのか
ここは誤解しやすいので丁寧に見ます。「数列を生成するLFSR」は一般に複数あります。極端な話、$L = 2t$ の長いLFSRなら、$2t$ 個の数列はいくらでも再現できてしまいます（初期値をそのまま並べればよい）。ではなぜ最短を選ぶのが正しいのでしょうか。
理由は、真の誤りパターンが最短解を与え、かつ $\nu \le t$ なら最短解が一意だからです。
実際に誤りが $\nu$ 個起きていれば、真の $\Lambda(x)$ は長さ $\nu$ のLFSRを与えます。もし別の $\Lambda'(x)$（長さ $L’ < \nu$）も同じ数列を生成できたとすると、$\Lambda'$ に対応する誤りパターンは $L'$ 個以下の誤りということになります。ところがRS符号の最小距離は $d_{\min} = 2t+1$ なので、重み $\nu + L' \le 2t$ の2つの誤りパターンが同じシンドロームを与えることはありません（差が符号語になるには重み $2t+1$ 以上が必要）。よって $\nu \le t$ の範囲では最短LFSRは真の誤りパターンに対応するものだけです。
逆にいうと、$\nu > t$ のときはこの保証が崩れます。$2t$ 個しかシンドロームがないので、長さ $t$ を超えるLFSRは同定できず、BMは「$L \le t$ の範囲で最善の嘘」を返します。後の実験でこの崩壊を実測します。
さて、最短LFSRを求めればよいことは分かりました。しかし「最短」を素直に探そうとすると、$L = 1$ から順に試して連立方程式を解く（$O(t^4)$ 相当）ことになり効率が悪い。BMアルゴリズムは、これをシンドロームを1個ずつ食べながら $O(t^2)$ でやってのけます。次節でその仕組みに入ります。
Berlekamp-Masseyアルゴリズム
発想: 増分的に設計図を直す
BMの基本方針は「とりあえず今までの分は合っているLFSRを持っておき、次の1個で食い違ったら最小限だけ直す」というものです。
イメージとしては、時系列予測モデルのオンライン学習に似ています。今のモデルで次の値を予測してみて、当たっていれば何もしない。外れたら、その外れ量に比例した補正を加える。BMがふつうのオンライン学習と違うのは、補正が完璧である点です。1回の補正で「過去すべて＋今回の1個」が厳密に合うようになります。近似ではなく代数的な厳密解を逐次的に構成していきます。
記号を用意します。$n = 0, 1, \dots, 2t-1$ をステップ番号とし（$S_{n+1}$ を処理するステップ）、次を保持します。
| 記号 | 意味 |
|---|---|
| $\Lambda^{(n)}(x)$ | ステップ $n$ 終了時点の接続多項式（$S_1 \dots S_{n+1}$ を生成できる） |
| $L$ | 現在のLFSRの長さ |
| $B(x)$ | 最後に長さを更新した直前の接続多項式（バックアップ） |
| $b$ | そのときのディスクレパンシー |
| $m$ | 最後に長さを更新してから経過したステップ数 |
初期値は $\Lambda^{(-1)} = 1$、$L = 0$、$B = 1$、$b = 1$、$m = 1$ です。長さ0のLFSRは「全部ゼロの数列」しか生成できないので、これが出発点になります。
ディスクレパンシー — 「どれだけ食い違ったか」
ステップ $n$ で、現在のLFSR $\Lambda^{(n-1)}$ に $S_{n+1}$ を予測させます。予測値は $-\sum_{i=1}^{L}\Lambda_i S_{n+1-i}$ ですから、実際の値との食い違いは
$$ \begin{equation} d_n = S_{n+1} + \sum_{i=1}^{L} \Lambda_i^{(n-1)} S_{n+1-i} \end{equation} $$
で測れます。これをディスクレパンシー（discrepancy、食い違い）と呼びます。$d_n = 0$ なら現在のLFSRは $S_{n+1}$ まで正しく生成できており、何もする必要がありません。
見方を変えると、$d_n$ は「$\Lambda^{(n-1)}(x) S(x)$ の $x^n$ の係数」です。鍵方程式(5)を思い出すと、$x^{\nu}$ 以降の係数はすべて0であるべきでした。$d_n$ はその「あるべきゼロ」からのズレを測っています。BMアルゴリズムとは、$\Lambda(x)S(x)$ の高次の係数を左から順にゼロに潰していく手続きだとも言えます。
更新則 — なぜ $x^m B(x)$ を足すのか
$d_n \ne 0$ のとき、$\Lambda$ を修正します。BMの更新則は次の形です。
$$ \begin{equation} \Lambda^{(n)}(x) = \Lambda^{(n-1)}(x) – \frac{d_n}{b} \, x^m \, B(x) \end{equation} $$
初見だと「なぜ $x^m B(x)$ なのか」がまったく分かりません。ここを丁寧に解きほぐします。
この図が更新則のすべてを説明しています。一番上の帯は $B(x)S(x)$ の係数列で、$m$ ステップ前の位置 $x^{n-m}$ にちょうど $b$ という値が立っています（それより左はゼロ、つまり $B$ はそこまで正しく生成できていた）。$x^m$ を掛けると、この帯全体が $m$ 個右にスライドし、$b$ が現在の位置 $x^n$ に到着します（中段）。下段は今の $\Lambda^{(n-1)}$ が出している食い違い $d_n$ です。$b$ を $d_n/b$ 倍して引けば $x^n$ の係数だけが正確にゼロになり、しかもシフトによって左側は依然ゼロのままなので過去の成果が一切壊れない — これが「1回の補正で厳密に合う」からくりです。
鍵は、$B(x)$ が「ちょうど $b$ という食い違いを出す装置」として保存されていることです。$B(x)$ は、いまから $m$ ステップ前（ステップ $n – m$）にディスクレパンシー $b \ne 0$ を出した接続多項式でした。式で書けば
$$ [x^{n-m}]\, B(x) S(x) = b $$
です。ここで $B(x)$ に $x^m$ を掛けると何が起きるでしょうか。多項式に $x^m$ を掛けることは、係数列を $m$ 個右にずらすことです。したがって
$$ [x^{n}]\, x^m B(x) S(x) = [x^{n-m}]\, B(x)S(x) = b $$
$m$ ステップ前に起きた食い違い $b$ が、ちょうど今のステップ $n$ に移動してくるのです。これが $x^m$ の役割です。
あとはスケールを合わせるだけです。$b$ を $d_n$ にしたいので、$d_n / b$ 倍します。すると
$$ [x^n]\left( \frac{d_n}{b} x^m B(x) S(x) \right) = \frac{d_n}{b} \cdot b = d_n $$
となり、これを $\Lambda^{(n-1)}$ から引けば
$$ [x^n]\, \Lambda^{(n)}(x)S(x) = d_n – d_n = 0 $$
で食い違いが消えます。まさに狙いどおりです。
さらに嬉しいことに、この補正は過去を壊しません。$x^m B(x)$ の係数は $x^m$ より低い次数がすべてゼロで、しかも $B(x)S(x)$ は $x^{n-m}$ より低い次数の係数がすべてゼロだったので（$B$ はそこまで正しく生成できていた）、$x^m B(x)S(x)$ は $x^{n}$ より低い次数の係数がすべてゼロです。つまり補正項は $x^n$ の係数だけをピンポイントで動かし、$x^0 \dots x^{n-1}$ には一切触れません。これが「1回の補正で厳密に合う」からくりです。
なお $\Lambda_0 = 1$（定数項）も保たれます。補正項の定数項は $m \ge 1$ よりゼロだからです。
長さ更新則 — いつレジスタを伸ばすか
もうひとつの要は、$L$ をいつ増やすかです。BMの規則は驚くほど簡潔です。
$$ \begin{equation} 2L \le n \ \text{ならば} \quad L \leftarrow n + 1 – L, \quad B \leftarrow \Lambda^{(n-1)}, \quad b \leftarrow d_n, \quad m \leftarrow 1 \end{equation} $$
$2L > n$ ならば $L$ は据え置きで、$m \leftarrow m+1$ とするだけです。
条件 $2L \le n$ の意味は「現在のLFSRが短すぎる」ということです。長さ $L$ のLFSRで $n+1$ 個の数列を説明しようとしていて、$n+1 > 2L$ なら情報が足りていません。このとき新しい長さは $n+1-L$ に跳ね上がります。次節でこれが最短であることを示します。
$B \leftarrow \Lambda^{(n-1)}$ という代入も重要です。いま失敗した $\Lambda^{(n-1)}$ そのものを、次回のための「ズレ発生装置」として保存するのです。失敗を捨てずに部品として再利用する、という設計思想が効いています。
アルゴリズム全体
以上をまとめると、BMアルゴリズムは次の擬似コードになります。
- 初期化: $\Lambda \leftarrow 1$, $B \leftarrow 1$, $L \leftarrow 0$, $b \leftarrow 1$, $m \leftarrow 1$
- $n = 0, 1, \dots, 2t-1$ について繰り返す:
- $d \leftarrow S_{n+1} + \sum_{i=1}^{L} \Lambda_i S_{n+1-i}$
- $d = 0$ なら $m \leftarrow m+1$ として次へ（何もしない）
- $2L \le n$ なら: $T \leftarrow \Lambda$; $\Lambda \leftarrow \Lambda – (d/b)x^m B$; $L \leftarrow n+1-L$; $B \leftarrow T$; $b \leftarrow d$; $m \leftarrow 1$
- そうでなければ: $\Lambda \leftarrow \Lambda – (d/b)x^m B$; $m \leftarrow m+1$
- $\Lambda(x)$ と $L$ を返す
ループは $2t$ 回、各回の計算は最大 $O(t)$ 回の有限体乗算なので、全体で $O(t^2)$ です。$t = 8$ なら16回のループ、乗算は高々200回程度。ハードウェアなら数十クロックで終わります。
「なぜこの更新則で最短性が保たれるのか」がまだ残っています。次節でMasseyの証明を追います。
最短性の証明 — Masseyの補題
何を示したいのか
$L_n$ を「$S_1, \dots, S_n$ を生成する最短LFSRの長さ」と定義します。示したいのは次の2点です。
(A) 下界: 長さ $L_n$ のLFSRが $S_{n+1}$ で失敗するなら、$L_{n+1} \ge n + 1 – L_n$ (B) 達成: BMの更新則が作る $\Lambda^{(n)}$ は長さ $n+1-L_n$ で、実際に $S_1, \dots, S_{n+1}$ を生成する
(B) は前節ですでに示しました。補正は過去を壊さず、$x^n$ の係数だけをゼロにするからです。長さについては、$\deg(x^m B) = m + \deg B$ が新しい長さの上界を与え、帰納法で $n+1-L$ に収まることが確認できます。難しいのは (A) の下界です。
下界の証明
背理法で示します。長さ $L = L_n$ のLFSR $\Lambda$ が $S_1, \dots, S_n$ を生成し、$S_{n+1}$ で失敗したとします。ここで、$S_1, \dots, S_{n+1}$ を生成する長さ $L’$ のLFSR $\Lambda’$ が存在し、$L’ \le n – L$ だと仮定します（示したい $L’ \ge n+1-L$ の否定）。
$\Lambda$ は $S_1, \dots, S_n$ を生成するので、$L+1 \le j \le n$ の範囲で
$$ S_j = -\sum_{i=1}^{L} \Lambda_i S_{j-i} $$
が成り立ちます。一方 $\Lambda’$ は $S_1, \dots, S_{n+1}$ を生成するので、$L’+1 \le j \le n+1$ の範囲で
$$ S_j = -\sum_{i=1}^{L’} \Lambda’_i S_{j-i} $$
が成り立ちます。
さて、$\Lambda$ が失敗した量 $d = S_{n+1} + \sum_{i=1}^{L}\Lambda_i S_{n+1-i} \ne 0$ を、$\Lambda’$ の漸化式を使って書き換えます。$S_{n+1}$ とその内側に現れる各 $S_{n+1-i}$ に、$\Lambda’$ の漸化式を代入していきます。
$$ d = S_{n+1} + \sum_{i=1}^{L}\Lambda_i S_{n+1-i} = \sum_{i=0}^{L} \Lambda_i S_{n+1-i} \quad (\Lambda_0 = 1) $$
各項に $\Lambda’$ の漸化式（$S_{n+1-i} = -\sum_{k=1}^{L’}\Lambda’_k S_{n+1-i-k}$、これは $n+1-i \ge L’+1$ すなわち $i \le n – L’$ で有効。$L’ \le n – L$ より $i \le L$ の全範囲で有効）を代入すると
$$ d = -\sum_{i=0}^{L}\sum_{k=1}^{L’} \Lambda_i \Lambda’_k S_{n+1-i-k} $$
となります。ここで和の順序を入れ替え、今度は内側の $\sum_i \Lambda_i S_{n+1-k-i}$ に $\Lambda$ の漸化式を適用します。
$$ d = -\sum_{k=1}^{L’} \Lambda’_k \left( \sum_{i=0}^{L} \Lambda_i S_{(n+1-k)-i} \right) $$
括弧の中は「$\Lambda$ に $S_{n+1-k}$ を予測させたときの食い違い」です。$1 \le k \le L’$ なので $n+1-k$ は $n+1-L’ \ge L+1$ 以上 $n$ 以下の範囲にあり、この範囲では $\Lambda$ は正しく生成できていました。つまり括弧の中はすべてゼロです。したがって $d = 0$ となり、$d \ne 0$ という仮定に矛盾します。
よって $L’ \ge n+1-L$ が示されました。つまり、失敗が起きた瞬間に、LFSRの長さは最低でも $n+1-L$ まで伸びざるを得ない。BMの更新則はこの下界をぴったり達成するので、常に最短解を保持し続けます。
直感的な読み方
この証明の心は「$\Lambda$ と $\Lambda’$ の両方が正しい領域が重なりすぎると、食い違い $d$ がゼロに追い込まれる」という点にあります。短いLFSRは「主張が強い」ぶん、過去のデータと矛盾しやすい。$L’ \le n-L$ という強い主張を置くと、$\Lambda$ の失敗そのものが打ち消されてしまうのです。
言い換えれば、BMは各ステップで「これ以上短くはできない」というギリギリの線を保ち続けます。貪欲に見えて、実は大域最適という性質は、アルゴリズム設計としても美しいところです。
理屈が固まったので、次は実際の数値で手を動かしてみましょう。
具体例 — GF(16) 上の RS(15,11) を手計算する
舞台設定
小さな体で全部の数値を追いかけます。GF(16) を原始多項式 $p(x) = x^4 + x + 1$ で構成します。$\alpha$ を $p(\alpha)=0$ の根とすると、$\alpha^4 = \alpha + 1$ という関係から冪表が作れます（16進表記、ビットが $\alpha^3\alpha^2\alpha^1\alpha^0$ の係数）。
| 冪 | 値 | 冪 | 値 | 冪 | 値 |
|---|---|---|---|---|---|
| $\alpha^0$ | 1 | $\alpha^5$ | 6 | $\alpha^{10}$ | 7 |
| $\alpha^1$ | 2 | $\alpha^6$ | C | $\alpha^{11}$ | E |
| $\alpha^2$ | 4 | $\alpha^7$ | B | $\alpha^{12}$ | F |
| $\alpha^3$ | 8 | $\alpha^8$ | 5 | $\alpha^{13}$ | D |
| $\alpha^4$ | 3 | $\alpha^9$ | A | $\alpha^{14}$ | 9 |
RS(15,11) は $n = 15$、$k = 11$、$2t = 4$、$t = 2$ です。生成多項式は
$$ g(x) = (x-\alpha)(x-\alpha^2)(x-\alpha^3)(x-\alpha^4) = x^4 + \alpha^{13}x^3 + \alpha^6 x^2 + \alpha^3 x + \alpha^{10} $$
となります（実際に計算すると係数は 16進で $1, \text{D}, \text{C}, 8, 7$）。
情報多項式を $u(x) = x^{10}$（つまり最上位シンボルだけが1、あとは0）として組織符号化すると、符号語は16進表記で
$$ c = \texttt{1 0 0 0 0 0 0 0 0 0 0 6 8 E 5} $$
となります（先頭11個が情報部、末尾4個がパリティ）。この符号語のシンドロームは当然すべて0です。
誤りを2個入れる
$x^{10}$ の位置に $\alpha^2$、$x^3$ の位置に $\alpha^5$ の誤りを加えます。つまり
$$ e(x) = \alpha^2 x^{10} + \alpha^5 x^3 $$
真の誤り位置数は $X_1 = \alpha^{10}$、$X_2 = \alpha^{3}$、誤り値は $Y_1 = \alpha^2$、$Y_2 = \alpha^5$ です。シンドロームを計算すると
$$ S_1 = \alpha^9, \quad S_2 = \alpha^8, \quad S_3 = \alpha^{13}, \quad S_4 = \alpha^7 $$
が得られます。検算しておきましょう。$S_1 = Y_1 X_1 + Y_2 X_2 = \alpha^2 \alpha^{10} + \alpha^5 \alpha^3 = \alpha^{12} + \alpha^8 = \text{F} \oplus 5 = \text{A} = \alpha^9$。確かに合っています。
BMアルゴリズムのトレース
初期状態は $\Lambda = 1$, $L = 0$, $B = 1$, $b = 1$, $m = 1$ です。
まず全体像を表で確認しておきます。列を左から追うと、各ステップで「$d_n$ を計算する → $2L \le n$ かどうかで長さを伸ばすか決める → $\Lambda$ を更新する」という同じ手順が4回繰り返されているだけだと分かります。注目すべきは4行目の $L$ で、$0 \to 1 \to 1 \to 2 \to 2$ と、長さ更新が起きた $n=0$ と $n=2$ でのみ増えています。1回の長さ更新に2個のシンドロームを消費するというリズムがすでにここに現れています。以下、各行を手計算で追いかけます。
$n = 0$（$S_1$ を処理）: $L = 0$ なので和は空、$d = S_1 = \alpha^9$。$d \ne 0$ かつ $2L = 0 \le 0$ なので長さ更新に入ります。$T = 1$、$\Lambda \leftarrow 1 – (\alpha^9/1)x^1 \cdot 1 = 1 + \alpha^9 x$。$L \leftarrow 0+1-0 = 1$、$B \leftarrow 1$、$b \leftarrow \alpha^9$、$m \leftarrow 1$。
$n = 1$（$S_2$ を処理）: $d = S_2 + \Lambda_1 S_1 = \alpha^8 + \alpha^9 \cdot \alpha^9 = \alpha^8 + \alpha^{18} = \alpha^8 + \alpha^3 = 5 \oplus 8 = \text{D} = \alpha^{13}$。$d \ne 0$ ですが $2L = 2 > 1 = n$ なので長さは据え置き、係数のみ補正します。$\Lambda \leftarrow (1 + \alpha^9 x) – (\alpha^{13}/\alpha^9)x \cdot 1 = 1 + (\alpha^9 + \alpha^4)x = 1 + (\text{A} \oplus 3)x = 1 + \alpha^{14}x$。$m \leftarrow 2$。
$n = 2$（$S_3$ を処理）: $d = S_3 + \Lambda_1 S_2 = \alpha^{13} + \alpha^{14}\alpha^8 = \alpha^{13} + \alpha^{22} = \alpha^{13} + \alpha^{7} = \text{D} \oplus \text{B} = 6 = \alpha^5$。$d \ne 0$、$2L = 2 \le 2 = n$ なので長さ更新です。$T = 1 + \alpha^{14}x$、$\Lambda \leftarrow (1+\alpha^{14}x) – (\alpha^5/\alpha^9)x^2 \cdot 1 = 1 + \alpha^{14}x + \alpha^{11}x^2$。$L \leftarrow 2+1-1 = 2$、$B \leftarrow 1+\alpha^{14}x$、$b \leftarrow \alpha^5$、$m \leftarrow 1$。
$n = 3$（$S_4$ を処理）: $d = S_4 + \Lambda_1 S_3 + \Lambda_2 S_2 = \alpha^7 + \alpha^{14}\alpha^{13} + \alpha^{11}\alpha^8 = \alpha^{10}$。$2L = 4 > 3$ なので係数のみ補正。$\Lambda \leftarrow (1+\alpha^{14}x+\alpha^{11}x^2) – (\alpha^{10}/\alpha^5)x(1+\alpha^{14}x) = 1 + \alpha^{12}x + \alpha^{13}x^2$。
シンドロームを使い切りました。結果は
$$ \Lambda(x) = 1 + \alpha^{12} x + \alpha^{13} x^2, \qquad L = 2 $$
です。$L = 2$ が「誤りは2個」と教えてくれています。
正解かどうかを検算する
真の誤り位置数から $\Lambda$ を組み立てると
$$ (1 – \alpha^{10}x)(1 – \alpha^{3}x) = 1 + (\alpha^{10} + \alpha^{3})x + \alpha^{13}x^2 $$
$\alpha^{10} + \alpha^3 = 7 \oplus 8 = \text{F} = \alpha^{12}$ なので、$1 + \alpha^{12}x + \alpha^{13}x^2$。BMの出力と完全に一致しました。4回のループで、非線形連立方程式を解いたのと同じ結果に到達したわけです。
Chien探索とForney公式で仕上げる
$\Lambda(x)$ の根を探します。位置 $i$ に誤りがある $\iff \Lambda(\alpha^{-i}) = 0$ なので、$i = 0, 1, \dots, 14$ をすべて代入します。実際に計算すると $i = 3$ と $i = 10$ でゼロになり、真の誤り位置と一致します。
誤り値は $\Omega(x)$ から求めます。$\Omega(x) = S(x)\Lambda(x) \bmod x^4$ を計算すると
$$ \Omega(x) = \alpha^9 + \alpha^{14} x $$
形式微分は、標数2では偶数次の項が消えるので $\Lambda'(x) = \Lambda_1 = \alpha^{12}$（定数）です。Forney公式（次節で導出）を使うと
$$ Y_1 = \frac{\Omega(\alpha^{-10})}{\Lambda'(\alpha^{-10})} = \alpha^2, \qquad Y_2 = \frac{\Omega(\alpha^{-3})}{\Lambda'(\alpha^{-3})} = \alpha^5 $$
となり、注入した誤り値 $\alpha^2, \alpha^5$ をぴったり復元できました。この小例だけで、復号の全パイプラインが検算済みになります。
手計算で流れをつかんだところで、誤り値の公式そのものをきちんと導いておきましょう。
Chien探索とForney公式
Chien探索 — 根の全数探索を効率化する
$\Lambda(x)$ の根を求める問題は、一般には難しい方程式です。しかし有限体では候補が $\alpha^0, \alpha^{-1}, \dots, \alpha^{-(n-1)}$ の $n$ 個しかないので、全部代入するのが最も確実です。これをChien探索（Chien search）と呼びます。
素朴に $\Lambda(\alpha^{-i})$ を毎回ホーナー法で計算すると $O(nL)$ の乗算が必要ですが、Chienの工夫は「隣の候補との差分」を使うことです。$\Lambda(\alpha^{-i}) = \sum_{j=0}^{L}\Lambda_j \alpha^{-ij}$ の各項 $\Lambda_j \alpha^{-ij}$ を保持しておき、$i \to i+1$ のときに各項に $\alpha^{-j}$ を掛けるだけで次の値が得られます。$L+1$ 個の乗算と加算で1候補が処理でき、ハードウェアでは $L+1$ 個の定数乗算器を並列に置いて1クロック1候補を実現します。
根が $L$ 個ちょうど見つかれば復号続行、足りなければ「訂正能力を超えた」と判断して復号失敗を宣言します。この失敗検出は実用上とても大事で、誤ったデータを黙って上位層に渡すよりずっと安全です。
Forney公式の導出
誤り位置が分かったので、あとは誤り値です。素朴には $\nu$ 元連立一次方程式（Vandermonde系）を解けばよいのですが、Forneyは割り算1回で済む閉じた式を与えました。
$\Omega(x)$ の定義式(4)に $x = X_k^{-1}$ を代入します。$l \ne k$ の因子は生き残り、$l = k$ の項だけが $\prod$ から外れているので
$$ \Omega(X_k^{-1}) = Y_k X_k \prod_{l \ne k}\left( 1 – X_l X_k^{-1} \right) $$
です。一方、$\Lambda(x) = \prod_l (1 – X_l x)$ の形式微分は、積の微分則から
$$ \Lambda'(x) = \sum_{l} (-X_l) \prod_{j \ne l} (1 – X_j x) $$
となります。ここで $x = X_k^{-1}$ を代入すると、$l \ne k$ の項には因子 $(1 – X_k X_k^{-1}) = 0$ が含まれるのですべて消え、$l = k$ の項だけが残ります。
$$ \Lambda'(X_k^{-1}) = -X_k \prod_{j \ne k}\left(1 – X_j X_k^{-1}\right) $$
2つの式を割り算すると、共通の積 $\prod_{l\ne k}(1 – X_l X_k^{-1})$ と $X_k$ がきれいに約分されます。
$$ \frac{\Omega(X_k^{-1})}{\Lambda'(X_k^{-1})} = \frac{Y_k X_k}{-X_k} = -Y_k $$
したがって
$$ \begin{equation} Y_k = -\frac{\Omega(X_k^{-1})}{\Lambda'(X_k^{-1})} \end{equation} $$
が得られました。これがForney公式です。標数2の体では $-1 = 1$ なので符号は無視できます。
なお、生成多項式の根を $\alpha^b, \alpha^{b+1}, \dots$（$b \ne 1$）と取る流儀では、$Y_k = -X_k^{1-b}\,\Omega(X_k^{-1})/\Lambda'(X_k^{-1})$ という補正因子が付きます。本記事では $b = 1$ を採るので補正は不要です。実装で他人のコードと結果が合わないときは、この $b$ の食い違いが原因であることが非常に多いので覚えておいてください。
形式微分の計算も標数2では簡単です。$\frac{d}{dx}x^j = j x^{j-1}$ で、$j$ が偶数なら $j \equiv 0$ なので項が消え、奇数なら $j \equiv 1$ で係数がそのまま残ります。つまり奇数次の項だけを1つ下にずらすだけで $\Lambda’$ が求まります。
理論は完結しました。ここからはPythonで最後まで動かします。
Pythonでの実装
GF($2^8$) の四則演算を対数表で作る
有限体の乗算をまじめに多項式乗算＋剰余で行うと遅いので、実装では対数表を使います。GF($2^8$)の乗法群は位数255の巡回群なので、原始元 $\alpha = 2$（多項式 $x$ に対応）の冪で全非零元を尽くせます。すると $ab = \alpha^{\log a + \log b}$ という具合に、乗算が整数の足し算に化けます。
原始多項式には $p(x) = x^8 + x^4 + x^3 + x^2 + 1$（16進で 0x11D）を使います。これはRS符号の標準的な選択です。
import numpy as np
# GF(2^8) の対数表・逆対数表を作る (原始多項式 x^8+x^4+x^3+x^2+1 = 0x11D)
PRIM = 0x11D
EXP = [0] * 512 # EXP[i] = alpha^i (二重長にして加算の折返しを省く)
LOG = [0] * 256 # LOG[v] = i such that alpha^i = v
_x = 1
for _i in range(255):
EXP[_i] = _x
LOG[_x] = _i
_x <<= 1 # x 倍 = 1ビット左シフト
if _x & 0x100: # 8次が立ったら原始多項式で還元
_x ^= PRIM
for _i in range(255, 512):
EXP[_i] = EXP[_i - 255] # alpha^255 = 1 なので周期255
def gmul(a, b):
"""GF(2^8) の乗算"""
if a == 0 or b == 0:
return 0
return EXP[LOG[a] + LOG[b]]
def gdiv(a, b):
"""GF(2^8) の除算 (b != 0)"""
if a == 0:
return 0
return EXP[(LOG[a] - LOG[b]) % 255]
def ginv(a):
"""乗法逆元"""
return EXP[(255 - LOG[a]) % 255]
def gpow(a, n):
"""冪乗 (n は負でもよい)"""
if a == 0:
return 0
return EXP[(LOG[a] * n) % 255]
# 加算・減算は XOR そのもの
print("alpha^0..8 =", [gpow(2, i) for i in range(9)])
print("分配則の確認:", gmul(0x53, 0x1B ^ 0xCA) == (gmul(0x53, 0x1B) ^ gmul(0x53, 0xCA)))
print("逆元の確認 :", all(gmul(a, ginv(a)) == 1 for a in range(1, 256)))
出力は alpha^0..8 = [1, 2, 4, 8, 16, 32, 64, 128, 29]、分配則・逆元の確認はともに True になります。$\alpha^8 = 29 = \text{0x1D}$ となるのは、$x^8 = x^4+x^3+x^2+1$（＝0x1D）という還元がここで初めて効くからです。1から128までは単なる倍々ですが、8乗で初めて原始多項式が顔を出す様子が数値としてはっきり見えます。逆元がすべての非零元で存在することも確認でき、GF($2^8$)が確かに体になっていることが実測できました。
多項式演算とRS(255,239)の符号化
多項式は昇順のリスト（インデックス＝次数）で表します。降順にすると添字が $\text{len}-1-i$ になってバグの温床になるので、BMのように「$i$ 次の係数」を頻繁に参照するアルゴリズムでは昇順が安全です。
# --- 多項式演算 (係数は昇順: p[i] が x^i の係数) ---
def pmul(p, q):
r = [0] * (len(p) + len(q) - 1)
for i, a in enumerate(p):
if a:
for j, b in enumerate(q):
if b:
r[i + j] ^= gmul(a, b)
return r
def padd(p, q):
n = max(len(p), len(q))
r = [0] * n
for i, c in enumerate(p):
r[i] ^= c
for i, c in enumerate(q):
r[i] ^= c
return r
def pscale(p, s):
return [gmul(c, s) for c in p]
def pshift(p, m):
"""x^m を掛ける = 係数を m 個右にずらす"""
return [0] * m + list(p)
def peval(p, xv):
"""ホーナー法で p(xv) を評価"""
y = 0
for c in reversed(p):
y = gmul(y, xv) ^ c
return y
def ptrim(p):
"""最高次側の余分な0を落とす"""
while len(p) > 1 and p[-1] == 0:
p = p[:-1]
return p
これでGF($2^8$)係数の多項式を自由に扱えます。peval にホーナー法を使っているのは、$L+1$ 回の乗算で評価が済むからです。Chien探索で $255 \times 9$ 回呼ぶことになるので、ここの効率は無視できません。
続いてRS(255,239)を組みます。$n = 255$、$k = 239$、$2t = 16$、$t = 8$ です。
N, K = 255, 239
NSYM = N - K # 2t = 16
T = NSYM // 2 # t = 8
def rs_generator(nsym):
"""g(x) = prod_{i=1}^{nsym} (x - alpha^i) を昇順係数で返す"""
g = [1]
for i in range(1, nsym + 1):
g = pmul(g, [gpow(2, i), 1]) # (alpha^i + x)
return g
G = rs_generator(NSYM)
def rs_encode(msg):
"""組織符号化。msg は長さ K の降順シンボル列 -> 長さ N の降順符号語"""
rem = [0] * NSYM # 剰余レジスタ
for c in msg:
f = c ^ rem[NSYM - 1] # フィードバック値
for i in range(NSYM - 1, 0, -1):
rem[i] = rem[i - 1] ^ gmul(G[i], f)
rem[0] = gmul(G[0], f)
return list(msg) + [rem[NSYM - 1 - i] for i in range(NSYM)]
def syndromes(cw):
"""降順符号語からシンドローム S_1..S_2t を計算"""
poly = list(reversed(cw)) # 昇順に直す
return [peval(poly, gpow(2, j)) for j in range(1, NSYM + 1)]
rng = np.random.default_rng(0)
msg = [int(v) for v in rng.integers(0, 256, K)]
cw = rs_encode(msg)
print("生成多項式の次数:", len(G) - 1)
print("無誤り符号語のシンドロームが全ゼロ:", not any(syndromes(cw)))
生成多項式の次数は16、無誤り符号語のシンドロームは全ゼロと表示されます。これは「符号語は必ず $g(x)$ の倍数であり、$g$ の根 $\alpha^1, \dots, \alpha^{16}$ で消える」という設計が正しく実装できたことの確認です。符号化はシフトレジスタによる剰余計算で、239シンボルを1回通すだけの $O(nk’)$ 相当の軽い処理です。
Berlekamp-Masseyの実装
理論どおりに書き下すと、驚くほど短いコードになります。
def berlekamp_massey(S, trace=False):
"""シンドローム列 S から誤り位置多項式 Lambda(昇順) と長さ L を返す"""
L = 0
Lam = [1] # 現在の接続多項式
B = [1] # 最後に長さ更新した直前のバックアップ
b = 1 # そのときのディスクレパンシー
m = 1 # 長さ更新からの経過ステップ数
rows = []
for n in range(len(S)):
# (1) ディスクレパンシー d = S[n] + sum_i Lam_i S[n-i]
d = S[n]
for i in range(1, L + 1):
if i < len(Lam): # Lam の次数が L 未満のこともある
d ^= gmul(Lam[i], S[n - i])
if d == 0:
m += 1
act = "d=0 (更新なし)"
else:
# (2) 補正項 (d/b) x^m B(x) を引く (標数2なので XOR)
corr = pshift(pscale(B, gdiv(d, b)), m)
new_Lam = ptrim(padd(Lam, corr))
if 2 * L <= n:
# (3) 長さ更新: L <- n+1-L, バックアップを取り替える
Lam, B, b, L, m = new_Lam, Lam, d, n + 1 - L, 1
act = "長さ更新"
else:
Lam = new_Lam
m += 1
act = "係数のみ補正"
rows.append((n, d, L, list(Lam), act))
return (Lam, L, rows) if trace else (Lam, L)
注意すべきは if i < len(Lam) のガードです。BMでは $\deg\Lambda < L$ となる場面があり得るため、$i$ が $\Lambda$ の長さを超えたときにインデックスがはみ出さないようにしなければなりません。ここを省くと、Pythonでは負インデックスが末尾に回り込んでまったく間違った $\Lambda$ を返すという発見しにくいバグになります（筆者は実際にこれで誤り3個の復号が1%失敗しました）。
Chien探索とForney公式
def chien_search(Lam):
"""Lambda(alpha^{-i}) = 0 となる位置 i を全数探索"""
return [i for i in range(N) if peval(Lam, gpow(2, (-i) % 255)) == 0]
def forney(S, Lam, positions):
"""Forney公式で誤り値を求める (b=1 の規約)"""
Om = pmul(S, Lam)[:NSYM] # Omega = S*Lam mod x^{2t}
Lp = [Lam[i] if i % 2 == 1 else 0 for i in range(1, len(Lam))] # 形式微分
if not Lp:
Lp = [0]
vals = {}
for p in positions:
Xinv = gpow(2, (-p) % 255)
den = peval(Lp, Xinv)
if den == 0: # 分母0 = 復号不能
return None
vals[p] = gdiv(peval(Om, Xinv), den)
return vals
def rs_decode(recv):
"""RS(255,239) の復号。(復号語, L, 成功フラグ) を返す"""
S = syndromes(recv)
if not any(S):
return list(recv), 0, True # 誤りなし
Lam, L = berlekamp_massey(S)
pos = chien_search(Lam)
if len(pos) != len(Lam) - 1: # 根が足りない = 訂正能力超過
return list(recv), L, False
vals = forney(S, Lam, pos)
if vals is None:
return list(recv), L, False
out = list(recv)
for p, e in vals.items():
out[N - 1 - p] ^= e # 冪 p は降順配列の index N-1-p
if any(syndromes(out)): # 訂正後の再検査
return out, L, False
return out, L, True
chien_search の結果の個数と $\deg\Lambda$ を突き合わせている点、訂正後にもう一度シンドロームを検査している点が、実用的な復号器で必須の安全弁です。訂正能力を超えたときにでたらめなデータを出さず、はっきり「失敗」と申告できます。
動かしてみる
誤りを2個入れて、全パイプラインが通るか確認します。
recv = list(cw)
recv[10] ^= 0x2C # 降順index 10 → 冪 244 の位置
recv[100] ^= 0x7A # 降順index 100 → 冪 154 の位置
S = syndromes(recv)
print("S1..S8 =", [hex(s) for s in S[:8]])
Lam, L, rows = berlekamp_massey(S, trace=True)
print("\n--- BMのトレース (最初の6ステップ) ---")
for n, d, LL, lam, act in rows[:6]:
print("n=%d d=0x%02X L=%d Lam=%s [%s]" % (n, d, LL, [hex(c) for c in lam], act))
pos = chien_search(Lam)
vals = forney(S, Lam, pos)
print("\nL =", L, " 誤り位置(冪) =", sorted(pos))
print("誤り値 =", {p: hex(v) for p, v in vals.items()})
dec, L, ok = rs_decode(recv)
print("復号成功:", ok, " 元の符号語と一致:", dec == cw)
実行すると次の出力が得られます。
S1..S8 = ['0xff', '0x27', '0x97', '0x35', '0xcc', '0x9c', '0xfe', '0x27']
--- BMのトレース (最初の6ステップ) ---
n=0 d=0xFF L=1 Lam=['0x1', '0xff'] [長さ更新]
n=1 d=0xC5 L=1 Lam=['0x1', '0x1f'] [係数のみ補正]
n=2 d=0x0D L=2 Lam=['0x1', '0x1f', '0x95'] [長さ更新]
n=3 d=0x25 L=2 Lam=['0x1', '0xc3', '0x54'] [係数のみ補正]
n=4 d=0x00 L=2 Lam=['0x1', '0xc3', '0x54'] [d=0 (更新なし)]
n=5 d=0x00 L=2 Lam=['0x1', '0xc3', '0x54'] [d=0 (更新なし)]
L = 2 誤り位置(冪) = [154, 244]
誤り値 = {154: '0x7a', 244: '0x2c'}
復号成功: True 元の符号語と一致: True
このトレースをグラフにすると、BMの挙動が視覚的に確認できます。上段のディスクレパンシーは $n=0,1,2,3$ の4ステップだけ非ゼロ（0xFF, 0xC5, 0x0D, 0x25）で、$n=4$ 以降は12ステップ連続でぴったりゼロです。$2\nu=4$ 個のシンドロームを食べた時点で答えが確定し、残りの12個は「その答えが正しい」ことを追認しているだけなのが読み取れます。下段の $L$ は $n=0$ と $n=2$ の2回だけ跳ね上がって2に到達し、そこで止まります。誤りが2個なら長さ2のLFSRで足りるので、BMはそれ以上むやみに伸ばそうとしません。
このトレースには本記事の理論がそのまま現れています。まず $n=0,1,2,3$ の4ステップだけでBMは正解の $\Lambda(x) = 1 + \text{0xC3}\,x + \text{0x54}\,x^2$ に到達し、$n=4$ 以降はディスクレパンシーが0のまま何も起きません。誤りが2個なら $2\nu = 4$ ステップで決まり、残り12個のシンドロームは「答え合わせ」に使われているだけだと分かります。次に、長さ更新が $n=0$ と $n=2$ の偶数ステップで起きています。これは $2L \le n$ という条件から自然に出てくるリズムで、1回の長さ更新に2個のシンドロームを消費する（未知数が位置と値で2個増えるから）という直感と一致します。最後に、誤り位置は冪 244 と 154 として返り、降順インデックス 10（$= 254-244$）と 100（$=254-154$）に対応しています。注入した 0x2C と 0x7A も完全に復元できました。
訂正能力の崖を実測する
RS(255,239) の理論上の訂正能力は $t = 8$ です。誤り数を0個から12個まで振って、訂正成功率がどう変化するかを実測します。
import matplotlib, matplotlib.pyplot as plt
for cand in ["Hiragino Sans", "Yu Gothic", "Noto Sans CJK JP", "IPAexGothic", "Meiryo"]:
if any(cand == f.name for f in matplotlib.font_manager.fontManager.ttflist):
plt.rcParams["font.family"] = cand
break
plt.rcParams["axes.unicode_minus"] = False
rng = np.random.default_rng(42)
trials = 400
err_counts = list(range(0, 13))
rates = []
for ne in err_counts:
success = 0
for _ in range(trials):
m_ = [int(v) for v in rng.integers(0, 256, K)]
c_ = rs_encode(m_)
r_ = list(c_)
if ne:
for p in rng.choice(N, size=ne, replace=False):
r_[int(p)] ^= int(rng.integers(1, 256)) # 非零の誤りを注入
d_, L_, ok_ = rs_decode(r_)
if ok_ and d_ == c_:
success += 1
rates.append(success / trials)
print("誤り数 %2d : 訂正成功率 %.3f" % (ne, success / trials))
plt.figure(figsize=(9, 5))
plt.plot(err_counts, rates, "o-", lw=2, ms=7, color="#1f77b4", label="訂正成功率（実測）")
plt.axvline(8, color="crimson", ls="--", lw=2, label="訂正能力 $t=8$")
plt.fill_between([-0.3, 8], 0, 1.05, color="#1f77b4", alpha=0.08)
plt.text(3.8, 0.55, "確実に訂正できる領域", ha="center", fontsize=12, color="#1f77b4")
plt.text(10.5, 0.55, "復号失敗を検出", ha="center", fontsize=12, color="crimson")
plt.xlabel("符号語1個あたりの誤りシンボル数 $\\nu$")
plt.ylabel("訂正成功率")
plt.title("RS(255,239) の訂正能力の崖（400試行/点）")
plt.ylim(-0.05, 1.08)
plt.xlim(-0.3, 12.3)
plt.grid(alpha=0.3)
plt.legend(loc="center left")
plt.tight_layout()
plt.show()
実行結果は、誤り数0〜8個ではすべて 1.000、9個以上ではすべて 0.000 という、教科書的なまでにくっきりした崖になります。
グラフにすると崖の鋭さがよく分かります。$\nu = 0$ から $8$ までは400試行すべてが成功して線が水平に張り付き、$\nu = 9$ で一気に床まで落ちて、以降は完全にゼロのままです。機械学習の分類器の性能曲線のようになだらかに劣化するのではなく、代数的な保証がある範囲は100%、その外は0%という二値的な振る舞いになっています。境界の位置も理論値 $t = (n-k)/2 = 8$ とぴったり一致しており、実装が正しいことの傍証にもなっています。
「成功率0.000」の中身も分解して確認しておきましょう。この棒グラフは $\nu = 9, 10, 11, 12$ の各400試行を「正しく訂正」「誤訂正（黙って誤ったデータを出す）」「復号失敗を自己申告」の3つに分類したものです。結果は1600試行すべてがオレンジ、つまり復号器は一度も嘘をつかず、必ず ok=False を返したことになります。赤（誤訂正）の面積はゼロです。ストレージや通信の実務では「訂正できなかった」より「間違ったデータを正しいと言い張った」ほうが遥かに危険なので、この性質は非常に重要です。
この結果から3つのことが読み取れます。第一に、訂正は確率的ではなく決定的です。誤りが8個以下なら、位置がどこであれ値が何であれ100%訂正できます。これは「最小距離 $2t+1$ の球充填」という代数的保証がそのまま出た形で、機械学習の分類器のように「だいたい当たる」のとは性質がまったく違います。第二に、9個になると成功率が0.5とか0.1に落ちるのではなく、いきなり0になります。訂正能力の外は完全に外なのです。第三に、失敗した400×4試行のすべてで復号器は ok=False（＝失敗の自己申告）を返し、誤ったデータを黙って出力するケース（誤訂正）は1件も起きませんでした。RS符号ではシンドロームに余裕があるため、多くの場合は誤訂正ではなく検出可能な失敗になります。
崖の正体をLFSRの長さで見る
なぜ9個で完全に破綻するのか。BMが返す長さ $L$ の推移を追うと理由が一目で分かります。
fig, ax = plt.subplots(figsize=(9, 5))
rng = np.random.default_rng(3)
for ne, style in [(1, "o-"), (3, "s-"), (5, "^-"), (8, "d-"), (12, "x--")]:
m_ = [int(v) for v in rng.integers(0, 256, K)]
c_ = rs_encode(m_)
r_ = list(c_)
for p in rng.choice(N, size=ne, replace=False):
r_[int(p)] ^= int(rng.integers(1, 256))
_, _, rows = berlekamp_massey(syndromes(r_), trace=True)
ax.plot(range(1, len(rows) + 1), [rr[2] for rr in rows], style,
lw=1.8, ms=6, label="誤り数 $\\nu=%d$" % ne)
ax.axhline(8, color="crimson", ls=":", lw=2)
ax.text(1.2, 8.25, "$L$ の上限 $t=8$", color="crimson", fontsize=11)
ax.set_xlabel("処理したシンドロームの個数 $n+1$")
ax.set_ylabel("そのときの最短LFSR長 $L$")
ax.set_title("BMアルゴリズムにおけるLFSR長の成長")
ax.set_xticks(range(1, 17))
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()
実際の $L$ の軌跡は次のようになります。
| 誤り数 $\nu$ | $L$ の軌跡（$n+1 = 1 \dots 16$） |
|---|---|
| 1 | 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1 |
| 3 | 1,1,2,2,3,3,3,3,3,3,3,3,3,3,3,3 |
| 5 | 1,1,2,2,3,3,4,4,5,5,5,5,5,5,5,5 |
| 8 | 1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8 |
| 12 | 1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8 |
このグラフは崖の正体を鮮やかに説明してくれます。まず、どのケースでも $L$ は2シンドロームにつき1ずつ階段状に増え、$\nu$ に達したところでぴたりと止まります。これはBMが「必要最小限まで伸ばし、それ以上は伸ばさない」最短性を実際に守っている証拠です。そして決定的なのは、$\nu=12$ の軌跡が $\nu=8$ の軌跡と完全に同一だという点です。シンドロームが16個しかないため、BMは長さ8を超えるLFSRを同定する材料を持たず、$L=8$ で止まらざるを得ません。返される $\Lambda$ は「16個のシンドロームには合うが、真の誤りとは無関係な偽物」です。だからChien探索で根が8個そろわず（あるいは根の個数が次数と食い違い）、復号失敗として検出されます。$t$ を超えた瞬間に成功率が0になるのは、確率の問題ではなく情報が足りていないという構造的な理由なのです。
理論から実装、そして限界の実測まで一周しました。最後に、BM以外の解き方や、実運用で組み合わされる拡張を眺めて、この手法の立ち位置を確かめておきましょう。
他手法との比較と発展
拡張ユークリッド互除法との関係
鍵方程式 $\Lambda(x)S(x) \equiv \Omega(x) \pmod{x^{2t}}$ は、$x^{2t}$ と $S(x)$ に拡張ユークリッド互除法を適用し、剰余の次数が $t$ 未満になった時点で止めることでも解けます。この方法をSugiyamaらのアルゴリズムと呼びます。
$$ \Lambda(x) x^{2t} + \Omega(x) \cdot (\cdots) = \text{（互除法の途中の剰余）} $$
という形で、$\Lambda$ が互除法のBézout係数として自然に現れます。BMとユークリッド法は数学的に等価で、実際「BMの各ステップはユークリッド法の1回の除算に対応する」ことが証明されています。
同じ鍵方程式を解く2つの流儀を並べると、得意分野が正反対なのが分かります。BMは状態を5つの変数しか持たず、シンドロームを1個ずつ流し込むだけなので、ハードウェアで1クロック1シンドロームのパイプラインに落としやすい構造です。一方ユークリッド法は多項式の除算を繰り返す分だけ一時的なデータ量が増えますが、$\Lambda$ と $\Omega$ が同時に手に入るのでForney公式にそのまま渡せます。どちらが「正しい」ということはなく、実装先がFPGA/ASICならBM、読みやすさ重視のソフトウェアならユークリッド法、という住み分けになります。実装上の使い分けを整理すると次のとおりです。
- BM: 逐次的でメモリが小さく、ハードウェア実装（特にパイプライン化）に向く。ただし体の除算 $d/b$ が毎ステップ必要（inversionless BMという除算を消した変種もある）
- ユークリッド法: $\Lambda$ と $\Omega$ が同時に手に入るので、Forney用の $\Omega$ を別計算しなくてよい。ソフトウェアでは分かりやすい
現代の高速実装では、ユークリッド法を分割統治で加速した $O(t \log^2 t)$ のアルゴリズムも知られていますが、$t$ が数十程度の実用域では定数倍で $O(t^2)$ のBMが勝つことが多いです。
消失訂正との併用
RS符号のもうひとつの強みは、位置が既知の誤り（消失、erasure）を効率よく扱えることです。例えば連結符号の内符号が「このブロックは怪しい」と申告した場合、位置は分かっているので未知数が半分になります。訂正可能条件は
$$ 2\nu + \rho \le 2t $$
（$\nu$ は位置未知の誤り数、$\rho$ は消失数）に緩みます。消失位置から消失位置多項式 $\Gamma(x) = \prod(1 – Z_j x)$ を作り、修正シンドローム $\tilde{S}(x) = \Gamma(x)S(x) \bmod x^{2t}$ にBMを適用すればよく、BM本体は一切変えずに済みます。
誤訂正の確率について
$\nu > t$ のとき、まれに「別の符号語に誤って復号される」誤訂正が起きます。本記事の実験では400×4試行で1件も観測されませんでしたが、確率は0ではありません。$t$ 個以下の誤りで到達できる符号語の球の総体積が符号空間全体に占める割合で概算でき、RS(255,239) では
$$ P_{\text{miscorrect}} \approx \frac{1}{t!} \approx \frac{1}{40320} \approx 2.5 \times 10^{-5} $$
というオーダーの見積もりが知られています。実験の1600試行では期待件数が0.04件なので、0件だったのは自然です。ミッションクリティカルな用途ではCRCを外側に重ねて、この残余リスクをさらに叩き落とします（誤り検出とCRC 参照）。
軟判定復号への橋渡し
BMは硬判定（受信シンボルを1つの値に確定してから復号する）の代表格です。近年は受信信号の尤度をそのまま使う軟判定復号が主流で、RS符号でもKötter-Vardyアルゴリズムのような代数的軟判定法や、リスト復号（Guruswami-Sudan）が研究されています。これらは $t$ を超える誤りでも一定の条件下で復号できますが、計算量は跳ね上がります。BMの $O(t^2)$ という軽さと決定性は、いまも組み込み用途で圧倒的な強みです。
代替手法と拡張まで見渡したところで、本記事でたどった道筋を整理しておきましょう。
まとめ
本記事では、Berlekamp-Masseyアルゴリズムを理論と実装の両面から解説しました。
- 鍵方程式: シンドロームの冪和表現から等比級数の和を経由して $\Lambda(x)S(x) \equiv \Omega(x) \pmod{x^{2t}}$ を導きました。$\deg\Omega < \deg\Lambda \le t$ という次数条件が、この合同式の解を一意にしています
- LFSRへの読み替え: 鍵方程式の $x^{\nu}$ 以降の係数比較から漸化式 $S_n = -\sum_i \Lambda_i S_{n-i}$ が出ます。これはLFSRの動作そのもので、「誤り位置を探す」問題が「シンドローム列を生成する最短LFSRを設計する」問題になります
- ディスクレパンシーと更新則: $d_n$ は $\Lambda(x)S(x)$ の $x^n$ 係数、つまり「あるべきゼロからのズレ」です。バックアップ $B(x)$ は「$m$ ステップ前に $b$ というズレを出した装置」なので、$x^m$ でずらして $d/b$ 倍すればズレをぴったり打ち消せます。しかも過去の係数には触れません
- 最短性: 失敗した瞬間に $L’ \ge n+1-L$ という下界が成立することを、2つのLFSRの漸化式を交互に代入して矛盾を導く形で証明しました。BMの更新則 $L \leftarrow n+1-L$ はこの下界をぴったり達成します
- 完全な実装: GF($2^8$)の対数表、RS(255,239)の符号化、BM、Chien探索、Forney公式をスクラッチで組み、誤り2個の復号を最後まで動かしました。BMは4ステップで正解にたどり着き、残り12個のシンドロームでは $d=0$ が続きます
- 訂正能力の崖: 誤り0〜8個で成功率1.000、9個以上で0.000という決定的な崖を実測しました。原因は「16個のシンドロームからは長さ8までのLFSRしか同定できない」ことにあり、$\nu=12$ の $L$ 軌跡が $\nu=8$ と完全一致することからも確認できます
BMアルゴリズムの美しさは、問題の言い換えにあります。非線形連立方程式という手強い姿を、「数列を生成する最小の回路」という機械的な問題に翻訳した瞬間に、$O(t^2)$ の逐次アルゴリズムが姿を現しました。似た構造は情報理論や信号処理のあちこちに顔を出します。Levinson-Durbinによる線形予測係数の高速計算、Padé近似、部分分数展開 — いずれも「短い漸化式で長い列を説明する」という同じ骨格を持っています。
次のステップとして、以下の記事も参考にしてください。
- リードソロモン符号の理論 — バースト誤りに強い符号を数学から理解する — 本記事で復号した符号の設計側
- BCH符号 — ガロア体上で設計する強力な誤り訂正符号 — 二元BCH符号ではシンドロームの共役性から奇数番だけを扱えばよくなり、BMがさらに簡略化されます
- 巡回符号の理論 — 生成多項式とシフトレジスタ符号化を基礎から解説 — 符号化側のシフトレジスタ実装
- LDPC符号の理論 — スパースグラフ符号の原理と復号アルゴリズム — 代数的復号とは対照的な、確率的・反復的な復号の世界
- 畳み込み符号とビタビ復号の理論 — トレリスと動的計画法 — RS符号と連結して使われる相方
- 誤り検出とCRC — パリティからCRC-32まで、データの正しさを守る技術 — 訂正ではなく検出に特化した、同じ多項式代数の応用