PyMatchingは、量子誤り訂正（QEC）符号の復号のための高速なPython/C++ライブラリで、MWPM（Minimum Weight Perfect Matching、最小重み完全マッチング問題）デコーダを使用しています。 量子誤り訂正回路からのシンドローム測定値が与えられると、 MWPMデコーダは、誤りの発生は独立であり、かつグラフ的であるという仮定（各誤りは1つまたは2つの検出イベントを引き起こす）のもと、最も可能性の高い誤り箇所の集合を探索します。 MWPMデコーダは、表面符号を復号するための最も一般的なデコーダです。 また、 サブシステム符号、 ハニカム符号、 2次元双曲線符号 といった、他のさまざまな符号のデコードに使用できます。
バージョン2にはblossomアルゴリズムの新しい実装が含まれており、以前のバージョンよりも100-1000倍高速になりました。 PyMatching は、境界のあるなしに関わらず、任意の重み付きグラフを使用して設定することができ、 Craig GidneyのStimライブラリと組み合わせて、回路レベルのノイズがある場合のエラー訂正回路のシミュレーションとデコードを行うことができます。 sinter パッケージでは、Stim と PyMatching を組み合わせて、量子エラー訂正回路のモンテカルロ・サンプリングを高速・並列に実行します。
PyMatchingのドキュメントはpymatching.readthedocs.io をご覧ください。
stim, sinter, pymatchingによって回路レベルのノイズを含む誤り訂正符号の閾値を推定する方法については、stim getting started notebookを試してみてください。
バージョン2は、私がCraig Gidneyと共同で書いたblossomアルゴリズムの新しい実装を特徴としています。 私たちの新しい実装は、sparse blossom アルゴリズムと呼ばれ、QECに関連する復号問題を扱うためのblossomアルゴリズムの一般化として見ることができます。 blossomアルゴリズムを一般化して、QECに関連する復号化問題を扱えるようにしたものです。 我々は、検出グラフの検出イベント間の最小重みの経路を見つける問題を直接解きます。 これによって、元のblossomアルゴリズムを使用した際に発生する派生グラフのMWPMを見つけるためのコストのかかる全対全ダイクストラ探索を回避することができます。 新しいバージョンはまたexactになっています。以前のバージョンのPyMatchingとは異なり、近似は行われません。
我々の新しい実装は、以前のバージョンのPyMatchingと比較して100倍以上高速です。 NetworkXに比べて100,000倍の速度です(表面符号回路でベンチマークした場合)。回路ノイズが0.1%の場合、PyMatchingは 距離13までの表面符号回路のXおよびZ基底を、シングルコアでのシンドローム解析1ラウンドあたり1マイクロ秒未満でデコードすることができます（X基底の測定値のみを処理する場合は距離19まで、ただし大規模な状況ではX基底とZ基底の両方を処理する必要があります）。 さらに、実行時間はグラフのノード数にほぼ比例しています。
以下のグラフでは、回路レベルの脱分極ノイズを含む表面符号回路のデコードにおいて、PyMatching v2と旧バージョン(v0.7)、NetworkXを比較したものです。 すべてのデコーダはM1プロセッサのシングルコアで実行され、XとZの両方の基底の測定値を処理しました。 凡例にあるT=N^xの式（および破線でプロット）は、同じデータに対するフィッティングによって得られたものです。 ここで、Nは1ラウンドあたりの検出器（ノード）数、Tは1ラウンドあたりの復号化時間です。 データおよびstim回路、他のベンチマークについてはリポジトリのbenchmarksフォルダを参照してください。
Sparse blossomはAustin Fowlerの論文に記載されているものと概念的に似たアプローチとなっていますが、我々のアプローチは細部が色々と異なっています(これについては我々の次の論文で説明する予定です)。 また、最近fusion-blossomライブラリーをリリースしたYue Wuによる素晴らしい独自研究とも類似点があります。 違いの1つは、fusion-blossomは交互木の探索領域を、Union-Findでクラスタが成長するのと同じような方法で成長させるのに対し、我々のアプローチは時系列に沿って進行することです。 そして、グローバルな優先順位キューを使用して、交互木を成長させます。 Yueは近々論文も発表する予定ですので、そちらもご期待ください。
PyMatchingの最新版はPyPIからダウンロードしてインストールすることができます。 以下のコマンドでインストールできます。
pip install pymatching --upgrade
PyMatching は stim.DetectorErrorModel 、 networkx.Graph 、 retworkx.PyGraph といったチェック行列からマッチンググラフを読み込むことができます。
あるいは、pymatching.Matching.add_edge と pymatching.Matching.add_boundary_edge を使って個別にエッジを追加することができます。
PyMatchingはStimと組み合わせることができます。
sinter (v1.10.0、あるいはより最新のバージョンを使用してください)は、PyMatchingとStimを使って、量子エラー訂正回路の並列モンテカルロシミュレーションを実行するもので、一般に最も簡単で速い方法です。
しかし、このセクションでは、StimとPyMatchingを直接使用し、Python APIの使用方法について説明します。
Stimをインストールするには、pip install stim --upgrade を実行してください。
まず、stimの回路を生成します。ここでは、stimに付属する表面符号回路を使用します。
import numpy as np
import stim
import pymatching
circuit = stim.Circuit.generated("surface_code:rotated_memory_x",
distance=5,
rounds=5,
after_clifford_depolarization=0.005)次に、stim を用いて stim.DetectorErrorModel (DEM) を生成します。これは、実質的には
Tannerグラフであり、回路レベルのノイズモデルを記述します。
decompose_errors=True を設定することにより、stim はすべてのエラーメカニズムを edge-like エラーメカニズム（1つまたは2つの検出イベントを引き起こす）に分解します。
これにより、DEMはグラフ的になり、pymatchingで読み込むことができるようになります。
model = circuit.detector_error_model(decompose_errors=True)
matching = pymatching.Matching.from_detector_error_model(model)次に、回路から1000ショットをサンプリングします。各ショット（shotsの列）は、完全なシンドローム（検出器
の測定値）と、ノイズの多い回路をシミュレートして得られた論理的な観測値からなります。
sampler = circuit.compile_detector_sampler()
syndrome, actual_observables = sampler.sample(shots=1000, separate_observables=True)これでデコードができるようになりました！PyMatchingの論理的観測値の予測結果と、stimでサンプリングした実際の観測値を比較し、誤りの数を数え、論理エラー率を推定します。
num_errors = 0
for i in range(syndrome.shape[0]):
predicted_observables = matching.decode(syndrome[i, :])
num_errors += not np.array_equal(actual_observables[i, :], predicted_observables)
print(num_errors) # prints 8また、バイナリのパリティチェック行列、Tannerグラフの別の表現、から pymatching.Matching オブジェクトをロードすることもできます。
パリティチェック行列 H の各行はパリティチェックに対応し、各列はエラーメカニズムに対応します。
H の要素 H[i,j] は、パリティチェック i がエラーメカニズム j によって反転された場合に 1 となり、そうでない場合に 0 となります。
PyMatching で使用するためには、H のエラーメカニズムは graphlike （グラフ的）である必要があります。
これは、各列には1つまたは2つの1が含まれていなければならないことを意味します（列が1つの1を持つ場合、それは境界に接続されたハーフエッジを表します）。
PyMatchingに weights というnumpy配列を与えることで、グラフ内の各エッジに重みを与えることができます。
weights 配列の要素 weights[j] は、 H の列 j に対応するエッジの重みを設定します。
エラー要因が独立したものとして扱われる場合、通常、エッジ j の重みは log((1-p_j)/p_j) として設定されます。
ここで、 p_j はエッジ j に関連する誤り確率です。
この設定により、PyMatchingはシンドロームが与えられたときに、最も確率の高いエラーメカニズムの集合を見つけます。
PyMatchingが H と weights を用いて設定されている場合、バイナリシンドロームベクトル syndrome (長さ H.shape[0] のnumpy配列) のデコードは、 H@predictions % 2 == syndrome を満たし、かつ解の総重み predictions@weights を最小にするような、バイナリ predictions ベクトルで定義されるエラーのセットを見つけることに対応します。
量子誤り訂正では、どのエラーメカニズムが発生したかを正確に予測するのではなく、通常、「論理的に観測可能な」測定結果、すなわちエラーメカニズムのパリティ、を予測したいのです。
これらはバイナリ行列の observables で表現されます。チェック行列と同様に、論理観測値 i がエラーメカニズム j によって反転された場合、observables[i,j] は 1 になります。
例えば、シンドローム syndrome が、一連のエラー noise （長さ H.shape[1] のバイナリ配列）の結果、たとえば、syndrome = H@noise % 2 であったとします。
もし、 observables@noise % 2 == observables@predictions % 2 ならば、我々のデコードは成功です。
これをまとめると、距離5の繰り返し符号を次のようにデコードすることができます。
import numpy as np
from scipy.sparse import csc_matrix
import pymatching
H = csc_matrix([[1, 1, 0, 0, 0],
[0, 1, 1, 0, 0],
[0, 0, 1, 1, 0],
[0, 0, 0, 1, 1]])
weights = np.array([4, 3, 2, 3, 4]) # Set arbitrary weights for illustration
matching = pymatching.Matching(H, weights=weights)
prediction = matching.decode(np.array([0, 1, 0, 1]))
print(prediction) # prints: [0 0 1 1 0]
# Optionally, we can return the weight as well:
prediction, solution_weight = matching.decode(np.array([0, 1, 0, 1]), return_weight=True)
print(prediction) # prints: [0 0 1 1 0]
print(solution_weight) # prints: 5.0また、物理エラー率10%に対する論理エラー率を推定するために、次のようにサンプルします。
import numpy as np
from scipy.sparse import csc_matrix
import pymatching
H = csc_matrix([[1, 1, 0, 0, 0],
[0, 1, 1, 0, 0],
[0, 0, 1, 1, 0],
[0, 0, 0, 1, 1]])
observables = csc_matrix([[1, 0, 0, 0, 0]])
error_probability = 0.1
weights = np.ones(H.shape[1]) * np.log((1-error_probability)/error_probability)
matching = pymatching.Matching.from_check_matrix(H, weights=weights)
num_errors = 0
for i in range(1000):
noise = (np.random.random(H.shape[1]) < error_probability).astype(np.uint8)
syndrome = H@noise % 2
prediction = matching.decode(syndrome)
predicted_observables = observables@prediction % 2
actual_observables = observables@noise % 2
num_errors += not np.array_equal(predicted_observables, actual_observables)
print(num_errors) # prints 4pymatching.Matching オブジェクト構築の際に faults_matrix 引数を与えることで、PyMatching に直接論理的な観測値を予測させることもできます。これによってデコーダは
いくつか追加の最適化を行い、デコード処理を少し速くすることができます。次の例では、この方法を使用しており、上で示した例と等価な操作になっています。
import numpy as np
from scipy.sparse import csc_matrix
import pymatching
H = csc_matrix([[1, 1, 0, 0, 0],
[0, 1, 1, 0, 0],
[0, 0, 1, 1, 0],
[0, 0, 0, 1, 1]])
observables = csc_matrix([[1, 0, 0, 0, 0]])
error_probability = 0.1
weights = np.ones(H.shape[1]) * np.log((1-error_probability)/error_probability)
matching = pymatching.Matching.from_check_matrix(H, weights=weights, faults_matrix=observables)
num_errors = 0
for i in range(1000):
noise = (np.random.random(H.shape[1]) < error_probability).astype(np.uint8)
syndrome = H@noise % 2
predicted_observables = matching.decode(syndrome)
actual_observables = observables@noise % 2
num_errors += not np.array_equal(predicted_observables, actual_observables)
print(num_errors) # prints 6チェック行列を使う代わりに、Matching オブジェクトをMatching.add_edge
メソッドと
Matching.add_boundary_edge
メソッドを用いて構築することもできます。または NetworkX や retworkx のグラフからロードすることによっても構築できます。
PyMatching の使い方の詳細については ドキュメントを参照してください。
PyMatching version 2で使用されている新しい実装（sparse blossom）に関する論文は近日中に発表される予定です。それまでの間、以下の論文の引用をお願いします。
@misc{pymatchingv2,
author = {Higgott, Oscar and Gidney, Craig},
title = {PyMatching v2},
year = {2022},
publisher = {GitHub},
journal = {GitHub repository},
howpublished = {\url{https://github.com/oscarhiggott/PyMatching}}
}
注釈: 既存のPyMatchingの論文では、バージョン0.7以前のPyMatchingの実装について記述されています。
We are grateful to the Google Quantum AI team for supporting the development of PyMatching v2. Earlier versions of PyMatching were supported by Unitary Fund and EPSRC.