量子情報理論の基本：データ圧縮（２）

$$
\def\bra#1{\mathinner{\left\langle{#1}\right|}}
\def\ket#1{\mathinner{\left|{#1}\right\rangle}}
\def\braket#1#2{\mathinner{\left\langle{#1}\middle|#2\right\rangle}}
$$

はじめに

前回、古典情報理論における「雑音のないチャネルにおけるシャノンの符号化定理」について理解できたので、今回はその量子版である「雑音のない量子チャネルにおけるシューマッハの符号化定理」について勉強します。シャノンの定理は、「典型系列」という概念を使って証明できましたが、シューマッハの定理も、類似の概念を使って、ほぼ同じ論理展開で証明できてしまいます。その説明の後、量子計算シミュレータqlazyを使って、符号化をシミュレーションしてみます。

参考にさせていただいたのは、以下の文献です。

1. ニールセン、チャン「量子コンピュータと量子通信（３）」オーム社（2005年）
2. 石坂、小川、河内、木村、林「量子情報科学入門」共立出版（2012年）
3. 富田「量子情報工学」森北出版（2017年）

エンタングルメント忠実度

ちょっとその前に、後の議論で使う「エンタングルメント忠実度(entanglement fidelity)」についておさえておきます。これまでの記事の中で一度も触れていなかったので、その定義と性質について、軽く確認しておきます。

エンタングメント忠実度というのは、CPTPマップで系の状態がどれだけ保持されるかを測る指標です。ある量子系Aにおける密度演算子$\rho^{A}$がCPTPマップによって$\Gamma(\rho^{A})$になったとします。このとき、$\rho^{A}$の純粋化を$\rho^{AR} = \ket{\Psi}^{AR}\bra{\Psi}^{AR}$とすると、エンタングルメント忠実度$F_{e}^{2}$は、2つの密度演算子同士の忠実度$F$を使って、

F_{e}^{2}(\rho^{A}, \Gamma) = F(\rho^{AR}, (\Gamma \otimes I^{R})(\rho^{AR}))^{2} = \bra{\Psi}^{AR} (\Gamma \otimes I^{R})(\rho^{AR}) \ket{\Psi}^{AR}  \tag{1}

のように定義されます¹。これは、純粋化した状態同士の内積に他ならないので、その上限は1になります。

また、CPTPマップはクラウス演算子$\{ M_k \}$を使って、$\Gamma(\rho) = \sum_{k} M_{k} \rho M_{k}^{\dagger}$と表現できました。とすると、エンタングルメント忠実度は、

F_{e}^{2}(\rho^{A}, \Gamma) = \sum_{k} |Tr(\rho^{A} M_k)|^{2}  \tag{2}

のようにも書けます。これは、

\begin{align}
F_{e}^{2}(\rho^{A}, \Gamma) &= \bra{\Psi}^{AR} \sum_{k} (M_k \otimes I_k) \rho^{AR} (M_{k}^{\dagger} \otimes I_k) \ket{\Psi}^{AR} \\
&= \sum_{k} \bra{\Psi}^{AR} (M_k \otimes I_k) \ket{\Psi}^{AR} \bra{\Psi}^{AR} (M_{k}^{\dagger} \otimes I_k) \ket{\Psi}^{AR}  \\
&= \sum_{k} Tr(\rho^{AR} (M_k \otimes I_k)) Tr(\rho^{AR} (M_{k}^{\dagger} \otimes I_k)) \\
&= \sum_{k} Tr(\rho^{A} M_k) Tr(\rho^{A} M_{k}^{\dagger}) \\
&= \sum_{k} |Tr(\rho^{A} M_k)|^{2}  \tag{3}
\end{align}

となることから、わかります。

問題設定

さて、それでは本題です。

今回、対象とするのは「古典情報」ではなく「量子情報」です。つまり、「古典情報」を「量子状態」に乗せて（＝符号化して）伝送し、受信側で測定（＝復号化して）「古典情報」を得るのではなく、量子状態そのものを伝送します。その際、いかにその「量子状態」を圧縮できるか、ということを考えていきます。と言葉で言ってしまうと簡単に聞こえるかもしれませんが、実際に物理的にどういうことをやろうとしているのか、まずはじめに確認しておきます。

量子情報理論において、量子状態と言えば、密度演算子$\rho$です。密度演算子で記述されるのは、一般に混合状態です。混合状態というのは、純粋状態の確率的なアンサンブルか、あるいは、大きな純粋状態の部分系として現れるものでした。例えて言うなら、前者はアリスが純粋状態を任意に設定できる光子発射装置を持っていて、それを確率的に切り替えて発射する際に現れる状態ですし、後者は量子コンピュータの一部の量子ビットを取り出すことで現れる状態です（例えば、、です）。

この密度演算子で記述される状態を一つの単位として、次々に伝送することを考えます。古典情報理論の場合、英語のアルファベットを時系列順に伝送することを、確率変数の系列$X_1,X_2,\cdots$でモデル化して表しましたが、それと同様に、伝送すべき量子状態を密度演算子の系列$\rho_1,\rho_2,\cdots$としてモデル化します。アリスの例で言うと、アリスはまず$\rho_1$に相当する光子発射行為を行い、一息ついて、次に$\rho_2$に相当する光子発射行為を行い、ということを繰り返すイメージです。あるいは、実は複数人のアリス（！）がいて、アリス1が$\rho_1$を担当し、アリス2が$\rho_2$を担当して、すべてのアリスが一斉に光子を発射するイメージでも良いかもしれません。量子コンピュータの例でいうと、ある量子回路で得られた量子状態の特定の部分量子ビットを取り出し、それを$\rho_1$と見なして、量子通信路に放出し、一息ついて、同じことをやって$\rho_2$を放出し、ということを繰り返して、密度演算子系列を生成するようなイメージです。あるいは、先程のアリスが複数人いる場合のように、大規模な量子回路を使って、同時並行に$\rho_1,\rho_2,\cdots$の系列を作っても良いかもしれません²。

ここで、生成された量子状態は、何らかの伝送路を介して伝送されます。アリスの光子は光ファイバーを通してボブのもとに届けられるかもしれませんし、量子コンピュータから放出された量子状態は、量子インターネットを介して別の量子コンピュータに届けられるかもしれません³。いずれにしろ、物理系が介在する話なので、瞬時に無限次元の量子状態を伝送できると思わない方が良いです。各々の$\rho_i$は、$d$次元のヒルベルト空間上（従って、量子ビット数は$\log d$）で定義されている密度演算子です。それを$n$個まとめて送信するとしたとき、情報生成時に必要になる量子ビット数は当然$n \log d$です。伝送路に物理的な制限があるはずなので、何らかの状態変換を施して、なるべくこの量子ビット数を削減した形で伝送し、受信側では受け取った量子状態に何らかの変換を施して、もともとの量子状態をなるべく誤差なく復元したいわけです。

このとき、圧縮できる量子ビット数に限界があるでしょうか？あるとしたら、それはどのように定量化されるものでしょうか？今回、勉強する「シューマッハの符号化定理」が、それに対する答えになります。

典型部分空間

まず、「シャノンの定理」のときと同様に、情報源をシンプルにモデル化します。各々の$\rho_i$は同じ確率分布をもち（＝行列表現したときに同じ要素をもつ）、異なる$i$に属する$\rho_i$は互いに独立である（＝相関がない）とします。このような量子情報源のことを「i.i.d量子情報源」と呼びます。以降の議論は、これを前提にします。

さて、密度演算子の系列$\{ \rho_i \}$をひとつの合成系$\sigma$と見なします。すなわち、

\sigma = \rho_1 \otimes \rho_2 \cdots  \tag{4}

です。この$\rho_1,\rho_2,\cdots$は区別のための添字がつけられていますが、実体は同じ形をしています。それを$\rho$という記号で代表させて、スペクトル分解すると、

\rho = \sum_{x} p(x) \ket{x} \bra{x}  \tag{5}

と書けます。$\rho$が1量子状態だとすると、そのスペクトル分解は直交する2つの状態、例えば、$\ket{up}\bra{up}$と$\ket{down}\bra{down}$についての線形和になり、その係数p(x)は状態$\ket{x}$が生起する確率と解釈できます。つまり、状態が$\ket{up}$である確率はp(up)、状態が$\ket{down}$である確率はp(down)で両者の和は1になります。また、このときの量子エントロピー$S(\rho)$は、p(x)={ p(up), p(down)}に対する古典エントロピー$H(p(x))$に等しいです。とすると、無限個の$\rho$の系列の中から$n$個を取り出して、それが状態$\ket{x_1}\ket{x_2} \cdots \ket{x_n}$である確率は、古典論の場合と全く同じように考えることができます⁴。すなわち、非常に生起しやすい典型的なパターンというのがあって、その確率は$p(x_1,x_2,\cdots,x_n) = p(x_1)p(x_2) \cdots p(x_n) \approx 2^{-nS(p)}$で近似できます。

そして、その典型的なパターンが生起する確率は、以下のような$\epsilon$の幅をもったものとして定量的に特徴づけることができます。

2^{-n(S(\rho)+\epsilon)} \leq p(x_1)p(x_2) \cdots p(x_n) \leq 2^{-n(S(\rho)-\epsilon)}  \tag{6}

これは、また、

|\frac{1}{n} \log \frac{1}{p(x_1) \cdots p(x_n)} - S(\rho) | \leq \epsilon  \tag{7}

と書き換えることができます⁵。

ここで、この$\epsilon$で特徴づけられる確率で生起する典型的な系列パターンに対応した状態$\ket{x_i}\ket{x_2} \cdots \ket{x_n}$の集合を「$\epsilon$典型部分空間($\epsilon$-typical subspace)」と呼び、$T(n,\epsilon)$と書くことにします。

とすると、$\epsilon$典型部分空間への射影演算子$P(n,\epsilon)$というものも定義できて、

P(n,\epsilon) = \sum_{x \in T(n,\epsilon)} \ket{x_1}\bra{x_1} \otimes \ket{x_2}\bra{x_2} \otimes \cdots \otimes \ket{x_n}\bra{x_n}  \tag{8}

のように書けます（このあたりは古典論にはない概念になります）。

さて、ここまで準備ができたところで、次に、典型部分空間に関して成り立つ３つの定理について見ていきます。

典型部分空間の定理(1)

ニールセン、チャンの記載をそのまま引用します。

「$\epsilon>0$を固定する。任意の$\delta>0$と十分大きな$n$に対して

Tr(P(n,\epsilon)\rho^{\otimes n}) \geq 1-\delta  \tag{9}

」

古典情報理論の「典型系列の定理(1)」とほぼ同じ内容です。前回の記事で、$\epsilon - \delta$の独特な言い回しが意味するところについては説明したので、ここでは繰り返しません。違いは、$\epsilon$典型部分空間である確率が$Tr$で表現されています。これは大丈夫でしょうか。もとの密度演算子を射影操作によって$\epsilon$典型部分空間へマッピングしておき、それのトレースをとっているので、これは$\epsilon$典型部分空間である確率になります。

では、証明してみます。

【証明】

\begin{align}
Tr(P(n,\epsilon) \rho^{\otimes n}) &= Tr(\sum_{x \in T(n,\epsilon)} p(x_1) \ket{x_1}\bra{x_1} \otimes p(x_2) \ket{x_2}\bra{x_2} \otimes \cdots \otimes p(x_n) \ket{x_n}\bra{x_n}) \\
&= \sum_{x \in T(n,\epsilon)} p(x_1) \cdots p(x_n) Tr(\ket{x_1}\bra{x_1} \otimes \cdots \otimes \ket{x_n}\bra{x_n}) \\
&= \sum_{x \in T(n,\epsilon)} p(x_1) \cdots p(x_n)  \tag{10}
\end{align}

ということで、式(9)の左辺は古典情報理論における「典型系列である確率」に等しいです。よって、前回の記事で証明した典型系列の定理(1)から、直ちに証明できます。（証明終）

典型部分空間の定理(2)

ニールセン、チャンの記載を引用します。

「任意の固定した$\epsilon>0$と$\delta>0$に対して$T(n,\epsilon)$の次元$|T(n,\epsilon)|=Tr(P(n,\epsilon))$は次式を満たす。

(1-\delta) 2^{n(S(\rho)-\epsilon)} \leq |T(n,\epsilon)| \leq 2^{n(S(\rho)+\epsilon)}  \tag{11}

」

ここで、$T(n,\epsilon)$の「数」ではなく「次元」と言っていることにご注意ください。古典情報において、$|T(n,\epsilon)|$は$\epsilon$典型系列の「数」を表していましたが、量子情報では、$T(n,\epsilon)$は、$\epsilon$典型部分空間、すなわち全体のヒルベルト空間に対する部分空間なので、数ではなく次元と呼ぶのが適切です。この部分空間へは射影演算子$P(n,\epsilon)$を施すことで到達できます。で、射影演算子$P(n,\epsilon)$は、式(8)で定義されていました。その$\sum$の中身は各々が直交する軸への射影であり、それを典型的になっているxの分だけ足し合わせる格好になっています。従って、古典情報の典型系列の数に相当するものは、量子情報においては、状態xのバリエーションの数、つまり、典型部分空間の次元$|T(n,\epsilon)$|ということになります。また、典型部分空間の次元は、式(8)の両辺のトレースをとればわかる通り、$Tr(P(n,\epsilon))$に等しいです。この定理(2)の記述に関する注意ポイントは以上です。それでは、証明します。

【証明】

前回の記事で証明した典型系列の定理(2)の$H(X)$を$S(\rho)$に置き換えることで、直ちに証明できます。（証明終）

典型部分空間の定理(3)

ニールセン、チャンの記載を引用します。

「$S(n)$はたかだか$2^{nR}$次元である$H^{\otimes n}$の任意の部分空間への射影オペレータとする。ここで$R<S(\rho)$は固定されている。このとき任意の$\delta>0$と十分大きな$n$に対して、

Tr(S(n)\rho^{\otimes n}) \leq \delta  \tag{12}

」

前回記事の典型系列の定理(3)と同様に、この定理の主張は、どんな小さな$\delta$をとったとしても、十分大きな$n$をとれば、式(12)が成り立つようにできる、ということです。要するに、$n \to \infty$の極限で$Tr(S(n)\rho^{\otimes n}) \to 0$となるということを言っています。では、証明します。

【証明】

$S(n)$を使った射影のトレース（つまり、確率の和）を、典型部分空間と非典型部分空間に関するトレースに分けます。

Tr(S(n)\rho^{\otimes n}) = Tr(S(n) \rho^{\otimes n} P(n,\epsilon)) + Tr(S(n) \rho^{\otimes n} (I-P(n,\epsilon)))  \tag{13}

ここで、$P(n,\epsilon)$と$\rho^{\otimes n}$は可換、すなわち、

\begin{align}
\rho^{\otimes n} P(n,\epsilon) &= \sum_{x \in T(n,\epsilon)} p(x_1) \cdots p(x_n) \ket{x_1}\bra{x_1} \cdots \ket{x_n}\bra{x_n} \\
&= P(n,\epsilon) \rho^{\otimes n}  \tag{14}
\end{align}

で、かつ、$(P(n,\epsilon))^2 = P(n,\epsilon)$なので（射影なので）、式(13)の右辺第1項は、

Tr(S(n) P(n,\epsilon) \rho^{\otimes n} P(n,\epsilon))  \tag{15}

となります。ここで、$P(n,\epsilon) \rho^{\otimes n} P(n,\epsilon)$の固有値は$p(x_1)p(x_2) \cdots p(x_n)$であり、式(6)よりその上限は$2^{-n(S(\rho)-\epsilon)}$です。式(15)のトレースの中の$S(n)$は$2^{nR}$個の固有値を取り出す操作に相当⁶し、トレースでそれらを足し合わせるので、結局、

Tr(S(n) P(n,\epsilon) \rho^{\otimes n} P(n,\epsilon)) \leq 2^{nR} 2^{-n(S(\rho)-\epsilon)} \approx 2^{-n(S(\rho)-R)}  \tag{16}

という不等式が成り立ち、$R<S(\rho)$という前提だったので、$n \to \infty$の極限で式(13)の右辺第1項は0に収束します。

次に、式(13)の右辺第2項について考えます。$S(n) \leq I$であり、$S(n)$と$\rho^{\otimes n}(I-P(n,\epsilon))$は両方とも正値なので、

0 \leq Tr(S(n)\rho^{\otimes n} (I-P(n,\epsilon))) \leq Tr(\rho^{\otimes n} (I-P(n,\epsilon)))  \tag{17}

が成り立ちます⁷。$n \to \infty$の極限で$I-P(n,\epsilon) \to 0$なので、$Tr(S(n)\rho^{\otimes n} (I-P(n,\epsilon))) \to \infty 0$、つまり、式(13)の右辺第2項も0に収束します。

従って、式(13)は、$n \to \infty$の極限で0に収束します。（証明終）

雑音のないチャネルにおけるシューマッハの符号化定理

さて、それでは本日のメイン・イベントであるシューマッハの定理を証明する準備が整いました。再び、ニールセン、チャンの記載を引用します。

「$\{ H,\rho\}$はi.i.d.量子情報源であるとする。もし$R>S(\rho)$ならば情報源$\{ H,\rho\}$に対してレート$R$の信頼性の高い圧縮法が存在する。もし$R<S(\rho)$ならばレート$R$の如何なる圧縮法も信頼性が低い。」

$\{ H,\rho\}$というのは、あるヒルベルト空間$H$上で定義された状態（密度演算子）$\rho$の系列のことで、それをi.i.d.量子情報源として、まずイメージしましょうというのが前提です。そうしたときに、伝送レート$R$を量子エントロピー$S(\rho)$よりも少しだけ大きい値にしておけば、信頼性が高い圧縮法が存在する、ということを主張しています。それでは、証明してみます。

【証明】

まず、$R>S(\rho)$の場合です。

$S(\rho)+\epsilon \leq R$を満たすような$\epsilon>0$を考えます。典型部分空間の定理(2)より、

|T(n,\epsilon)| \leq 2^{n(S(\rho)+\epsilon)} \leq 2^{nR}  \tag{18}

が成り立ちます。ここで、$H_{c}^{n}$は任意の$2^{nR}$次元のヒルベルト空間で、$T(n,\epsilon)$を含むものとします⁸。

符号化は次のようにやります。互いに直交する射影演算子$P(n,\epsilon), I-P(n,\epsilon)$で量子情報源から生成された量子情報$\sigma \equiv \rho_1 \otimes \cdots \rho_n$を測定します。得られた結果を"0","1"と表すことにすると、”0"が得られた場合は、典型部分空間上の状態だったことになりますし、"1"だった場合は、そうではなく非典型部分空間上の状態だったという解釈になります。"0"だった場合、量子状態をそのまま送り、"1"だった場合は何か適当な状態に置き換えて送ります（例えば、$\ket{00 \cdots 0}$とか）。

この結果、もとの量子状態は、以下のように変換されることになります。

C^{n}(\sigma) \equiv P(n,\epsilon) \sigma P(n,\epsilon) + \sum_{i} A_i \sigma A_i^{\dagger} \tag{19}

ここで、$C^{n}$は、

C^{n}: H^{\otimes n} \rightarrow H_{c}^{n}  \tag{20}

のように写像する演算です。また、右辺第2項は、非典型部分空間だった場合、適当な状態に置き換えるという効果を記述しており、例えば、

A_i \equiv \ket{0}\bra{i}  \tag{21}

のようなものと思えば良いです（$\ket{i}$は典型部分空間の直交補空間に対する基底です）。

復号の演算$D^{n}$は、恒等演算とします⁹。

D^{n}(\sigma) = \sigma  \tag{22}

もとの状態が符号化と復号化の過程で、どれだけ変化したかをエンタングルメント忠実度$F_{e}^{2}$で評価します。エンタングルメント忠実度の最大値が1であることと、クラウス演算子による式(2)を使い、さらに、典型部分空間の定理(1)より、

\begin{align}
1 &\geq F_{e}^{2}(\rho^{\otimes n}, D^{n} \circ C^{n}) \\
&\geq |Tr(\rho^{\otimes n} P(n,\epsilon))|^{2} \\
&\geq |1-\delta|^{2} \geq 1-2\delta  \tag{23}
\end{align}

が成り立ちます。従って、どんなに小さい$\delta>0$をとっても、十分大きな$n$をもってくれば上式を成り立たせることができます。つまり、$n \to \infty$の極限で、エンタングルメント忠実度は1に収束します。つまり、レート$R$の信頼できる圧縮法が存在するということが言えます。

次に、$R<S(\rho)$の場合です。

符号化演算は、$H^{\otimes n}$から$2^{nR}$次元の部分空間への射影演算子$S(n)$による写像と考えて良いです。

符号化演算$C^{n}$に対するクラウス演算子を$\{ C_j \}$、復号化演算$D^{n}$に対するクラウス演算子を$\{ D_k \}$とおきます。すなわち、符号化は、

\begin{align}
&\rho_1 \otimes \cdots \otimes \rho_n \\
&\rightarrow \rho_{1}^{\prime} \otimes \cdots \otimes \rho_{n}^{\prime} \\
&= C_{1} \rho_1 C_{1}^{\dagger} \otimes \cdots \otimes C_{n} \rho_n C_{n}^{\dagger} = \sum_{j} C_{j} \rho^{\otimes n} C_{j}^{\dagger}  \tag{24}
\end{align}

で、復号化は、

\begin{align}
&\rho_{1}^{\prime} \otimes \cdots \otimes \rho_{n}^{\prime} \\
&\rightarrow \rho_{1}^{\prime\prime} \otimes \cdots \otimes \rho_{n}^{\prime\prime} \\
&= D_{1} \rho_{1}^{\prime} D_{1}^{\dagger} \otimes \cdots \otimes D_{n} \rho_{n}^{\prime} D_{n}^{\dagger} = \sum_{k} D_{k} \rho^{\prime \otimes n} D_{k}^{\dagger}  \tag{25}
\end{align}

です。とすると、符号化・復号化によるエンタングルメント忠実度は、

F_{e}^{2}(\rho^{\otimes n}, D^{n} \circ C^{n}) = \sum_{j,k} |Tr(D_{k} C_{j} \rho^{\otimes n})|^{2}  \tag{26}

となります。$C_j$は、$S(n)$によって射影される$2^{nR}$次元の部分空間と同じ部分空間への射影演算なので、

C_j = S(n) C_j  \tag{27}

です。また、$S(n)$で写像された$2^{nR}$次元のヒルベルト空間が、$D_k$で写像される部分空間への射影演算子を$S^{k}(n)$とすると、

S^{k}(n) D_{k} S(n) = D_{k} S(n)  \tag{28}

が成り立ちます。式(27)と式(28)より、

D_{k} C_{j} = D_{k} S(n) C_{j} = S^{k}(n) D_{k} S(n) C_{j} = S^{k}(n) D_{k} C_{j} \tag{29}

これを、式(26)に代入すると、

F_{e}^{2}(\rho^{\otimes n}, D^{n} \circ C^{n}) = \sum_{j,k} |Tr(D_{k} C_{j} \rho^{\otimes n} S^{k}(n))|^{2}  \tag{30}

となります。さらに、コーシー・シュワルツの不等式から、

F_{e}^{2}(\rho^{\otimes n}, D^{n} \circ C^{n}) \leq \sum_{j,k} Tr(D_{k} C_{j} \rho^{\otimes n} C_{j}^{\dagger} D_{k}^{\dagger}) Tr(S^{k}(n) \rho^{\otimes n}) \tag{31}

が成り立ちます¹⁰。

典型部分空間の定理(3)より、右辺の2番目のトレースは、任意の$\delta>0$に対して、

Tr(S^{k}(n) \rho^{\otimes n}) \leq \delta  \tag{32}

とできます。従って、式(31)は、

F_{e}^{2}(\rho^{\otimes n}, D^{n} \circ C^{n}) \leq \delta \sum_{j,k} Tr(D_{k} C_{j} \rho^{\otimes n} C_{j}^{\dagger} D_{k}^{\dagger}) = \delta \tag{33}

となります。ここで、最後の等号は、$C^{n},D^{n}$によって密度演算子のトレースは保存することから成り立ちます。式(33)より、$n \to \infty$の極限で、エンタングルメント忠実度は0に収束するので、結局、$R<S(\rho)$の場合、如何なる圧縮法も信頼性が低いということが証明されました。（証明終）

符号化シミュレーション

ニールセン、チャンの「コラム12.4」にシューマッハの符号化の例が出ているので、それを実装してみたいと思います。1量子ビット状態として、

\rho = \frac{1}{4}
\begin{pmatrix}
3 & 1 \\
1 & 1
\end{pmatrix}  \tag{34}

という密度演算子を考えます。これを4つまとめて伝送します。つまり、$\rho^{\otimes 4}$を一つの系列として伝送します。シューマッハの符号化定理によると、伝送レートは量子エントロピー$S(\rho)$より少し大きいものを考えておけば良かったです。式(34)の場合、計算すると$S(\rho)$は約0.6なので、伝送レートを$R=0.75$として、4量子ビットを3量子ビットにまで圧縮する想定にします。で、復号側で再び4量子ビット状態を復元するということをやってみたいと思います。

シューマッハの符号化に際し、典型部分空間に射影するのですが、まず、その前に、$\rho$を対角化する基底で考える必要があります。つまり、$\rho$をスペクトル分解します（式(5)のように）。

\rho = p \ket{\bar{0}} \bra{\bar{0}} +(1-p) \ket{\bar{1}} \bra{\bar{1}}  \tag{35}

となります。ここで、

\begin{align}
\ket{\bar{0}} &= \cos \frac{\pi}{8} \ket{0} + \sin \frac{\pi}{8} \ket{1}  \\
\ket{\bar{1}} &= -\sin \frac{\pi}{8} \ket{0} + \cos \frac{\pi}{8} \ket{1}  \tag{36}
\end{align}

p = \frac{1}{4} (3 + \tan \frac{\pi}{8}) = 0.85..  \tag{37}

です¹¹。

このような基底で考えると、$\rho$は、状態$\ket{\bar{0}}$となる確率が$p \approx 0.85$、$\ket{\bar{1}}$となる確率が$1-p \approx 0.15$となる純粋状態の混合となります。もとの$\rho$の式(34)は基底$\ket{0},\ket{1}$で表現されているので、シューマッハの符号化に際し、まずはじめにやらないといけないのは、このような基底変換に相当するユニタリ変換です。

式(36)を見れば簡単にわかる通り、そのユニタリ変換はY軸まわりに$-\pi/4$だけ回転させる回転ゲートです。それを4量子ビットの各々別個に適用すれば良いです。

そうすると、$\ket{\bar{0}}$が生起する確率は$p \approx 0.85$なので、この回転ゲートを適用した結果の状態

\sum_{x_1=\bar{0}}^{\bar{1}} \sum_{x_2=\bar{0}}^{\bar{1}} \sum_{x_3=\bar{0}}^{\bar{1}} \sum_{x_4=\bar{0}}^{\bar{1}} C(x_1,x_2,x_3,x_4) \ket{x_1 x_2 x_3 x_4}  \tag{38}

において、$\ket{x_1 x_2 x_3 x_4}$の中に登場する$\bar{0}$の数が多いほど、その確率が大きいと言えます¹²。

この$\bar{0}$の数が多い状態というのは典型部分空間を構成する状態なので、これを判別して、適切な量子状態を割り当てて伝送すれば良いのですが、この操作をユニタリ変換として実行しないといけないというのが、量子情報理論の難しさです。考え方は、$\bar{1}$の数が少ないものから順に並べ、それを辞書的に順序付けるような対応で変換させます。具体的に列挙すると、以下の通りです（いちいち数字の上にバーを付けるのが面倒なので、ここでは省略しました）。

\begin{align}
& \ket{0000} \rightarrow \ket{0000} \\
& \ket{0001} \rightarrow \ket{0001} \\
& \ket{0010} \rightarrow \ket{0010} \\
& \ket{0100} \rightarrow \ket{0011} \\
& \ket{1000} \rightarrow \ket{0100} \\
& \ket{0011} \rightarrow \ket{0101} \\
& \ket{0101} \rightarrow \ket{0110} \\
& \ket{0110} \rightarrow \ket{0111} \\ 
& \ket{1001} \rightarrow \ket{1000} \\
& \ket{1010} \rightarrow \ket{1001} \\
& \ket{1100} \rightarrow \ket{1010} \\
& \ket{0111} \rightarrow \ket{1011} \\
& \ket{1011} \rightarrow \ket{1100} \\
& \ket{1101} \rightarrow \ket{1101} \\
& \ket{1110} \rightarrow \ket{1110} \\
& \ket{1111} \rightarrow \ket{1111} \tag{39}
\end{align}

という変換を実行すれば良いです。どうでしょうか。$\bar{1}$の数が少ないほど、上の方に順序付けられていて、つまり確率的にありがちな状態たちが上の方に固まっていて、特定の量子ビット（4つあるうちの右側3つの量子ビット）だけで、多くの状態を表現できるように変換がなされています。これは全体として置換演算になっているので、ユニタリです。つまり、各量子ビットに基底変換を適用した後、この置換に相当するユニタリ変換を実行すれば良いです。

伝送レートは、3量子ビットとしたので、一番左の量子ビットを捨てて、右側3つの量子ビットのみを残して伝送します。

復号する際には、まず、捨ててしまった量子ビット分を$\ket{0}$で追加して、4量子ビット状態にして、式(39)の置換の逆変換を実行します。最後に、4量子ビットの各々に対して、基底変換の逆変換を実施します。これで、もとの状態がほぼ戻っているはずです。

実装

以上の流れを実際にやってみました。全体のPythonコードを示します。

【2021.9.5追記】qlazy最新版でのソースコードはここに置いてあります。

import numpy as np
from qlazypy import QState,DensOp

PERM = {'0000':'0000', '0001':'0001', '0010':'0010', '0100':'0011',
        '1000':'0100', '0011':'0101', '0101':'0110', '0110':'0111', 
        '1001':'1000', '1010':'1001', '1100':'1010', '0111':'1011',
        '1011':'1100', '1101':'1101', '1110':'1110', '1111':'1111'}

def perm_matrix(PERM):

    dim = len(PERM)
    perm_mat = np.zeros((dim,dim))
    iperm_mat = np.zeros((dim,dim))

    for k,v in PERM.items():
        col = int(k,2)
        row = int(v,2)
        perm_mat[row][col] = 1
        iperm_mat[col][row] = 1

    return (perm_mat, iperm_mat)

def create_densop():

    mat_1 = np.array([[0.75,0.25],[0.25,0.25]])
    de_1 = DensOp(matrix=mat_1)
    de_ini = de_1.composite(num=4)
    de_1.free()

    return de_ini

def coder(de_ini, id_all, id_comp, theta, perm_mat):

    de_tmp = de_ini.clone()
    [de_tmp.ry(i, phase=-theta) for i in id_all]
    de_tmp.apply(perm_mat)
    de_comp = de_tmp.partial(id_comp)  # 4-qubit -> 3-qubit state
    de_tmp.free()

    return de_comp

def decoder(de_comp, id_all, theta, iperm_mat):

    qs_0 = QState(1)
    de_0 = DensOp(qstate=[qs_0])
    de_fin = de_0.tenspro(de_comp)  # add 1-qubit (|0>)
    de_fin.apply(iperm_mat)
    [de_fin.ry(i, phase=theta) for i in id_all]

    qs_0.free()
    de_0.free()

    return de_fin

if __name__ == '__main__':

    # settings
    id_all = [0,1,2,3]
    id_comp = [1,2,3]
    theta = 0.25
    (perm_mat,iperm_mat) = perm_matrix(PERM)

    # initial state (4-qubit state)
    de_ini = create_densop()

    # coding (4-qubit -> 3-qubit state)
    de_comp = coder(de_ini, id_all, id_comp, theta, perm_mat)

    # decoding (3-qubit -> 4-qubit state)
    de_fin = decoder(de_comp, id_all, theta, iperm_mat)

    # fidelity
    print("fidelity = {:.6f}".format(de_fin.fidelity(de_ini)))

    de_ini.free()
    de_comp.free()
    de_fin.free()

何をやっているか、ざっくり説明します。main処理部を見てください。

# settings
id_all = [0,1,2,3]
id_comp = [1,2,3]
theta = 0.25
(perm_mat,iperm_mat) = perm_matrix(PERM)

で、全体の量子ビット番号のリスト(id_all)と圧縮対象とする量子ビット番号のリスト(id_comp)を設定し、基底変換のための角度(theta)を$0.25\pi$としています¹³。さらに、関数perm_matrixで、置換演算のための変換行列とその逆変換行列を作成しています。その際、一番上で定義している辞書データPERMを参照しています（式(39)に対応した辞書です）。処理内容は関数定義を見てください。この置換行列は、後で状態に対して適用するユニタリ変換を作るときに使用します。

# initial state (4-qubit state)
de_ini = create_densop()

で、初期状態（密度演算子）を作成します。関数create_densopの中では、式(34)に対応した1量子の密度演算子を4つ集めて積状態を作る処理を実施しています。

de_ini = de_1.composite(num=4)

という1行がその処理です。compositeメソッド(v0.0.31で追加)で引数4を与えることで4量子の積状態を返してくれます。

# coding (4-qubit -> 3-qubit state)
de_comp = coder(de_ini, id_all, id_comp, theta, perm_mat)

で、符号化します。4量子ビットの状態を3量子ビットに圧縮します。関数coderの中身を見てください。

def coder(de_ini, id_all, id_comp, theta, perm_mat):

    de_tmp = de_ini.clone()
    [de_tmp.ry(i, phase=-theta) for i in id_all]
    de_tmp.apply(perm_mat)
    de_comp = de_tmp.partial(id_comp)  # 4-qubit -> 3-qubit state
    de_tmp.free()

    return de_comp

1行目で、初期状態(de_ini)をcloneメソッドでテンポラリ(de_tmp)にコピーしておいて、これに対して処理を実行します。2行目で、回転ゲートで基底変換します(v0.0.30で、密度演算子に対するゲート演算処理を追加しました。量子ゲート$U$に対して、$U \rho U^{\dagger}$を実行します)。内包表記を使えば1行で済みます。3行目で、置換に相当するユニタリ変換を実行します。applyメソッドで置換行列を指定すれば、これも1行で済みます。本当は量子ゲートの組み合わせで実行したいところですが、端折りました、汗¹⁴。最後、4行目で、一番左側の量子ビット（プログラムの中では0番目の量子ビット）を捨てて（部分トレースをとって）、3量子ビット状態にします(de_comp)。これをリターンします。

main処理部に戻ります。

# decoding (3-qubit -> 4-qubit state)
de_fin = decoder(de_comp, id_all, theta, iperm_mat)

で、3量子ビット状態を入力し、復号化した4量子ビット状態を出力します。関数decoderの中身を見てください。

def decoder(de_comp, id_all, theta, iperm_mat):

    qs_0 = QState(1)
    de_0 = DensOp(qstate=[qs_0])
    de_fin = de_0.tenspro(de_comp)  # add 1-qubit (|0>)
    de_fin.apply(iperm_mat)
    [de_fin.ry(i, phase=theta) for i in id_all]

    qs_0.free()
    de_0.free()

    return de_fin

1行目、2行目で追加する量子状態$\ket{0}$に対する密度演算子を用意します(de_0)。3行目で、その量子状態(de_0)と圧縮した量子状態(de_comp)の積状態を作ります(de_fin)。4行目で、置換演算の逆を実行します。5行目で、各量子ビットに基底変換の逆を施します。最後に、その結果(de_fin)をリターンします。

main処理部に戻ります。

# fidelity
print("fidelity = {:.6f}".format(de_fin.fidelity(de_ini)))

もとの状態(de_ini)と復号結果の状態(de_fin)の忠実度を計算して、表示します。以上です。

動作確認

実行結果を示します。忠実度だけ表示しているので、とてもシンプルです。忠実度が1.0であれば、誤差なしの圧縮実現ということなのですが、

fidelity = 0.970213

となり、少し誤差が出ました。

式(39)を見れば、明らかですが、下の方に行けば行くほど、めったに生起しない状態になります。とは言え、有限な確率値を持っているはずです。その下半分をごっそり削除するような圧縮をやっているので、誤差が出るのは当然です。シューマッハの符号化定理が主張するのは、$n \to \infty$の極限で、この誤差が十分に小さくなるということです。なので、$n$の値を$40,400,4000,\cdots$と変えてみると、誤差が無限に小さくなる様子が確認できるのでしょう。が、シミュレータでそんな巨大な量子ビット数を試せないので、$n=4$の場合の結果をもって、シューマッハの符号化シミュレーションができた、ということにしておきます。

おわりに

今回、「雑音のない場合」について勉強したので、次は「雑音のある場合」と行きたいところですが、その前に「量子誤り訂正」について理解しておいた方が多分良いので、次回は、そっち方面を攻めてみたいと思います。いろいろ取り上げるべき話題がありそうなので、何回か費やすと思います。その後、「雑音のある場合」とか、その上で「量子暗号」かなー、というのが、とりあえずの予定です。

以上