$$
\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}}
$$
はじめに
量子状態の違いを測る指標として、前回の記事では「フィデリティ」を取り上げたので、今回はもう一つの指標である「トレース距離(trace distance)」を勉強してみます。その定義と性質について説明した後、量子計算シミュレータqlazyを使って、その重要な性質について、実際に計算して確認してみたいと思います。
参考にさせていただいたのは、以下の文献です。
トレース距離とは
定義
ある状態$\rho$と$\sigma$があったとき、トレース距離は、以下のように定義されます。
D(\rho,\sigma) = \frac{1}{2} ||\rho-\sigma|| = \frac{1}{2} Tr|\rho-\sigma| \tag{1}
$\rho-\sigma$は（正値ではない）エルミートなので、
D(\rho,\sigma) = \frac{1}{2} Tr \sum_{i} |\lambda_{i}| \ket{i} \bra{i} = \frac{1}{2} \sum_{i} |\lambda_{i}| \tag{2}
のように、$\rho-\sigma$の固有値の絶対値の和として計算することができます1。
また、射影演算子または正値演算子$P$を使って、以下のように表すこともできます。
D(\rho,\sigma) = \max_{P} Tr(P(\rho-\sigma)) \tag{3}
この表式から、トレース距離というのは、状態$\rho$と状態$\sigma$があったときに、考えられるすべての測定（射影またはPOVM）を実施したときの確率の差の最大値を表している、という物理的な解釈が成り立ちます。この式(3)は、次節の性質を議論する際に何度も使う大事な式なので、証明したいのですが、その前に、任意のエルミート演算子$A$に対して、以下が成り立つことを証明しておく必要がありますので、まずこっちを先にやっつけておきます。
Tr(A_{+}) = \max_{0 \leq P \leq I} (PA) \tag{4}
ここで、$P$は その値が1以下になる $I$以下である正値演算子であり、
\begin{align}
A &= \sum_{i} a_{i} \ket{i} \bra{i} \\
A_{+} &= \sum_{i (a_i > 0)} a_{i} \ket{i} \bra{i} \\
A_{-} &= \sum_{i (a_i \leq 0)} |a_{i}| \ket{i} \bra{i} \tag{5}
\end{align}
です。$A_{+}, A_{-}$をエルミート演算子$A$の正部分(positive part)、負部分(negative part)と呼びます。
【証明】
\begin{align}
Tr(PA_{-}) &= Tr(P \sum_{i(a_i \leq 0)} |a_i| \ket{i} \bra{i}) \\
&= \sum_{i(a_i \leq 0),j} \bra{j} |a_i| P \ket{i} \braket{i}{j} \\
&= \sum_{i(a_i \leq 0)} |a_i| \bra{i} P \ket{i} \geq 0 \tag{6}
\end{align}
$A=A_{+}-A_{-}$なので、
Tr(PA) = Tr(PA_{+}) - Tr(PA_{-}) \leq Tr(PA_{+}) \leq Tr(A_{+}) \tag{7}
が成り立ちます。最後の不等号は、$I-P \geq 0$より、$Tr((I-P)A_{+})\geq 0$となることを使いました。これより、
Tr(A_{+}) = \max_{0 \leq P \leq I} Tr(PA) \tag{8}
が成り立ちます。ここで、$PA=A_{+}$を満たす$P$のときに最大値をとります。（証明終）
それでは、式(3)の証明です。
【証明】
エルミート演算子$A$に対して、
\begin{align}
A &= A_{+} - A_{-} \\
|A| &= A_{+} + A_{-} \tag{9}
\end{align}
なので、
\begin{align}
Tr(A_{+}) &= \frac{1}{2} (Tr|A| + Tr(A)) \\
Tr(A_{-}) &= \frac{1}{2} (Tr|A| - Tr(A)) \tag{10}
\end{align}
です。式(10)および$\rho-\sigma$のトレースが0であることを使って、
\begin{align}
D(\rho,\sigma) &= \frac{1}{2} Tr|\rho-\sigma| \\
&= Tr(\rho-\sigma)_{+} - \frac{1}{2} Tr(\rho-\sigma) \\
&= Tr(\rho-\sigma)_{+} \tag{11}
\end{align}
これに、式(8)を代入して、
D(\rho,\sigma) = \max_{0 \leq P \leq I} Tr(P(\rho-\sigma)) \tag{12}
となります。（証明終）
性質
トレース距離には、以下に示す性質があります2。
- (1)対称性： $D(\rho,\sigma) = D(\sigma,\rho)$
- (2)非負性(最大値は1)： $0 \leq D(\rho,\sigma) \leq 1$
- (3)三角不等式： $D(\rho,\sigma) \leq D(\rho,\tau) + D(\tau,\sigma)$
- (4)収縮性： $D(\rho,\sigma) \geq D(\Gamma(\rho),\Gamma(\sigma))$
- (5)強凸性(strongly convexity)： $D(\sum_{i} p_{i} \rho_{i}, \sum_{i} q_{i} \sigma_{i}) \leq D(p_{i},q_{i}) + \sum_{i} D(\rho_{i},\sigma_{i})$
この(1)(2)(3)は「距離」が満たすべき３つの条件です3。(4)は、物理過程でトレース距離は必ず減少する、つまり、周囲の環境との相互作用によって状態は次第に識別できなくなっていくということを表しています。
それでは、順に確認します。
(1)対称性
D(\rho,\sigma) = D(\sigma,\rho) \tag{13}
が成り立ちます。これは、定義から明らかですね。
(2)非負性(最大値は1)
0 \leq D(\rho,\sigma) \leq 1 \tag{14}
が成り立ちます。特に、最大・最小になるのは、以下のような場合です。
\begin{align}
\rho = \sigma &\Leftrightarrow D(\rho,\sigma) = 0 \\
\rho \sigma = 0 &\Leftrightarrow D(\rho,\sigma) = 1 \tag{15}
\end{align}
【証明】
まず、左側の不等号については、定義より明らかです（式(2)を見てください）。また、等号が成り立つのは、$\rho=\sigma$の場合のみです。
次に、右側の不等号についてです。$\rho,\sigma$を以下のように純粋化します。
\begin{align}
\rho &\rightarrow \ket{\phi_{\rho}} \\
\sigma &\rightarrow \ket{\phi_{\sigma}} \tag{16}
\end{align}
一般に、純粋化した状態間のトレース距離は、純粋化前のトレース距離よりも大きいので4、
\begin{align}
D(\rho,\sigma) &\leq \frac{1}{2} || \ket{\phi_{\rho}} \bra{\phi_{\rho}} - \ket{\phi_{\sigma}} \bra{\phi_{\sigma}} || \\
&= \frac{1}{2} Tr | \ket{\phi_{\rho}} \bra{\phi_{\rho}} - \ket{\phi_{\sigma}} \bra{\phi_{\sigma}} | \tag{17}
\end{align}
が言えます。上限が1であることを示すため、右辺のトレースを具体的に計算したいわけですが、そのためには何らかの正規直交系がほしくなります。そこでグラム・シュミットの直交化法を使います。ランクは高々２なので、最初の２つの直交系が得られれば良いです。というわけで、やってみます。
\begin{align}
\ket{u} &= \ket{\phi_{\rho}} \\
\ket{v} &= \frac{\ket{\phi_{\sigma}} - \braket{\phi_{\rho}}{\phi_{\sigma}} \ket{\phi_{\rho}}}{\sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}}} \tag{18}
\end{align}
したがって、
\begin{align}
\ket{\phi_{\rho}} &= \ket{u} \\
\ket{\phi_{\sigma}} &= \braket{\phi_{\rho}}{\phi_{\sigma}} \ket{u} + \sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}} \ket{v} \tag{19}
\end{align}
これより、式(17)のトレースの中身は（地道に計算すると）、
\begin{align}
&\ket{\phi_{\rho}} \bra{\phi_{\rho}} - \ket{\phi_{\sigma}} \bra{\phi_{\sigma}} \\
&= (1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}) \ket{u} \bra{u} - \sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}} \braket{\phi_{\rho}}{\phi_{\sigma}} \ket{u} \bra{v} - \sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}} \braket{\phi_{\sigma}}{\phi_{\rho}} \ket{v} \bra{u} - (1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}) \ket{v} \bra{v} \tag{20}
\end{align}
となり、$\alpha=\braket{\phi_{\rho}}{\phi_{\sigma}}$とおくと、
\begin{align}
&\ket{\phi_{\rho}} \bra{\phi_{\rho}} - \ket{\phi_{\sigma}} \bra{\phi_{\sigma}} \\
&= (1-|\alpha|^{2}) \ket{u}\bra{v} - \alpha \sqrt{1-|\alpha|^{2}} \ket{u}\bra{v} - \alpha^{*} \sqrt{1-|\alpha|^{2}} \ket{v}\bra{u} - (1-|\alpha|^{2}) \ket{v}\bra{v} \\
&=
\begin{pmatrix}
1-|\alpha|^{2} & - \alpha \sqrt{1-|\alpha|^{2}} \\
- \alpha^{*} \sqrt{1-|\alpha|^{2}} & - (1-|\alpha|^{2})
\end{pmatrix} \tag{21}
\end{align}
です。この行列を対角化し（＝固有値問題を解き）、得られた固有値の絶対値の和を計算すれば良いです。というわけで、やっみます。固有方程式は、求めたい固有値を$\lambda$として、
\lambda^{2} - (1-|\alpha|^{2}) = 0 \tag{22}
となります。これを解くと２つの固有値は、
\begin{align}
&\lambda_1 = \sqrt{1-|\alpha|^{2}} = \sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}} \\
&\lambda_2 = -\sqrt{1-|\alpha|^{2}} = -\sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}} \tag{23}
\end{align}
となり、結局、式(17)は、
D(\rho,\sigma) \leq \sqrt{1-|\braket{\phi_{\rho}}{\phi_{\sigma}}|^{2}} \leq 1 \tag{24}
となります。ここで、等号が成り立つのは、$\braket{\phi_{\rho}}{\phi_{\sigma}}=0$、すなわち状態が直交しているとき、すなわち、$\rho \sigma = 0$のときのみです。（証明終）
(3)三角不等式
D(\rho,\sigma) \leq D(\rho,\tau) + D(\tau,\sigma) \tag{25}
が成り立ちます。
【証明】
式(3)より、
D(\rho,\sigma) = Tr(P(\rho-\sigma)) \tag{26}
を満たす$P$が存在します。適当な$\tau$を使って、
D(\rho,\sigma) = Tr(P(\rho-\sigma)) = Tr(P(\rho-\tau)) + Tr(P(\tau-\sigma)) \tag{27}
とできます。再び、式(3)を使って、
D(\rho,\sigma) \leq D(\rho,\tau) + D(\tau,\sigma) \tag{25}
が成り立ちます。（証明終）
(4)収縮性
D(\rho,\sigma) \geq D(\Gamma(\rho),\Gamma(\sigma)) \tag{28}
が成り立ちます。
【証明】
注目系Aにおける状態$\rho,\sigma$に対して、参照系Rの状態$\tau_{R}$を使って、以下のように純粋化します。
\rho_{AR} = \rho \otimes \tau_{R}, \space \sigma_{AR} = \sigma \otimes \tau_{R} \tag{29}
このように同じ$\tau_{R}$を使って純粋化すると、トレース距離は不変であることに注意しておきます。すなわち、
\begin{align}
D(\rho_{AR}, \sigma_{AR}) &= \frac{1}{2} ||\rho_{AR}-\sigma_{AR}|| = \frac{1}{2} ||(\rho-\sigma) \otimes \tau_{R}|| \\
&= \frac{1}{2} ||\rho-\sigma|| \cdot ||\tau_{R}|| = \frac{1}{2} ||\rho-\sigma|| \\
&= D(\rho,\sigma) \tag{30}
\end{align}
です。また、$\rho,\sigma$に対するCPTPマップは、任意のユニタリ演算子$U$を使って、
\begin{align}
\Gamma(\rho) &= Tr_{R} (U \rho_{AR} U^{\dagger}) \\
\Gamma(\sigma) &= Tr_{R} (U \sigma_{AR} U^{\dagger}) \tag{31}
\end{align}
と書けます。そうすると、
\begin{align}
||\rho-\sigma|| &= ||\rho_{AR}-\sigma_{AR}|| = ||U(\rho_{AR}-\sigma_{AR})U^{\dagger}|| \\
&= \max_{V} |Tr(U(\rho_{AR}-\sigma_{AR})U^{\dagger}V)| \\
&\geq \max_{V_A} |Tr(U(\rho_{AR}-\sigma_{AR})U^{\dagger}(V_A \otimes I_R))| \\
&= \max_{V_A} |Tr((Tr_{R}(U\rho_{AR}U^{\dagger}) - Tr_{R}(U\sigma_{AR}U^{\dagger}))V_A)| \\
&= \max_{V_A} |Tr((\Gamma(\rho)-\Gamma(\sigma))V_A)| \\
&= ||\Gamma(\rho)-\Gamma(\sigma)|| \tag{32}
\end{align}
となります。
ここで、線形演算子$A$とすべてのユニタリ演算子$V$に対して、
||A|| = \max_{V} |Tr(AV)| \tag{33}
が成り立つことを使いました。結局、式(32)より、
D(\rho,\sigma) \geq D(\Gamma(\rho),\Gamma(\sigma)) \tag{28}
が成り立つことがわかります。（証明終）
(5)強凸性
D(\sum_{i} p_{i} \rho_{i}, \sum_{i} q_{i} \sigma_{i}) \leq D(p_{i},q_{i}) + \sum_{i} D(\rho_{i},\sigma_{i}) \tag{34}
が成り立ちます。
【証明】
\begin{align}
&D(\sum_{i} p_{i} \rho_{i}, \sum_{i} q_{i} \sigma_{i}) \\
&= \frac{1}{2} Tr|\sum_{i} p_{i} \rho_{i} - \sum_{i} q_{i} \sigma_{i}| \\
&= \max_{0 \leq P \leq I} Tr(P(\sum_{i} p_{i} \rho_{i} - \sum_{i} q_{i} \sigma_{i})) \tag{35}
\end{align}
ここで、$P$を適当に選ぶことで、以下のようにできます。
\begin{align}
&D(\sum_{i} p_{i} \rho_{i}, \sum_{i} q_{i} \sigma_{i}) \\
&= Tr(P(\sum_{i} p_{i} \rho_{i} - \sum_{i} q_{i} \sigma_{i})) \\
&= Tr(P \sum_{i} p_{i} \rho_{i}) - Tr(P \sum_{i} q_{i} \sigma_{i}) \\
&= \sum_{i} p_{i} Tr(P\rho_{i}) - \sum_{i} q_{i} Tr(P\sigma_{i}) \\
&= \sum_{i} p_{i} Tr(P(\rho_{i}-\sigma_{i})) + \sum_{i} (p_i - q_i) Tr(P\sigma_{i}) \\
&\leq \sum_{i} p_{i} D(\rho_{i},\sigma_{i}) + \sum_{i} (p_i - q_i) Tr(P\sigma_{i}) \tag{36}
\end{align}
ここで、最終行の第２項を古典確率分布の分布間距離に関連づけるため、$\rho,\sigma$が可換だった場合のトレース距離について考えます。可換なので、以下のように同時対角化が可能です。
\begin{align}
\rho &= \sum_{i} p_{i} \ket{i} \bra{i} \\
\sigma &= \sum_{i} q_{i} \ket{i} \bra{i} \tag{37}
\end{align}
この条件の元で、トレース距離を計算すると、
\begin{align}
D(\rho,\sigma) &= \frac{1}{2} Tr|\sum_{i} (p_i - q_i) \ket{i} \bra{i}| \\
&= \frac{1}{2} \sum_{i} |p_i - q_i| \equiv D(p_i,q_i) \\
&= \max_{P} Tr(P(\rho-\sigma)) \\
&= \max_{P} Tr(P(\sum_{i} p_{i} \ket{i} \bra{i} - \sum_{i} q_{i} \ket{i} \bra{i})) \\
&= \max_{P} \sum_{i} (p_i - q_i) Tr (P \ket{i} \bra{i}) \\
&= \max_{P} \sum_{i} (p_i - q_i) Tr (P \sigma_{i}) \\
&\geq \sum_{i} (p_i -q_i) Tr(P\sigma_{i}) \tag{38}
\end{align}
が導けます。ここで、$D(p_i,q_i)$は古典確率分布の分布間距離です。式(38)を式(36)に代入すると、
D(\sum_{i} p_{i} \rho_{i}, \sum_{i} q_{i} \sigma_{i}) \leq D(p_i,q_i) + \sum_{i} p_{i} D(\rho_{i},\sigma_{i}) \tag{34}
が成り立つことがわかります。（証明終）
ちなみに、$p_i = q_i$の場合、
D(\sum_{i} p_{i} \rho_{i}, \sum_{i} p_{i} \sigma_{i}) \leq \sum_{i} p_{i} D(\rho_{i},\sigma_{i}) \tag{39}
が成り立ち、これを「結合凸性(jointly convexity)」と呼んでいます。
シミュレータで確認
それでは、上で示したトレース距離の性質のうち、４番目の「収縮性」に注目し、それが本当なのかどうかをシミュレータで確認してみます。具体的には、２つの密度演算子をランダムに作成して、ランダムに作った量子チャネル(参照系＋環境系を追加して純粋化した状態に対してランダムなユニタリ変換を実施し最後にトレースアウトする、というやり方で定義しました)を通した結果、トレース距離が確かに収縮（＝減少）ということを見てみます。
全体のPythonコードは以下です。
【2021.9.5追記】qlazy最新版でのソースコードはここに置いてあります。
import random
import math
import numpy as np
from scipy.stats import unitary_group
from qlazypy import QState, DensOp
def random_densop(qnum_tar,qnum_ref,qnum_env):
dim_pur = 2**(qnum_tar+qnum_ref)
vec_pur = np.array([0.0]*dim_pur)
vec_pur[0] = 1.0
mat_pur = unitary_group.rvs(dim_pur)
vec_pur = np.dot(mat_pur, vec_pur)
dim_env = 2**qnum_env
vec_env = np.array([0.0]*dim_env)
vec_env[0] = 1.0
vec_whole = np.kron(vec_pur,vec_env)
qs = QState(vector=vec_whole)
de = DensOp(qstate=[qs],prob=[1.0])
qs.free()
return de
def random_unitary(qnum):
dim = 2**qnum
mat = unitary_group.rvs(dim)
return mat
if __name__ == '__main__':
# settings
qnum_tar = 2 # system A : target system
qnum_ref = 2 # system R : reference system
qnum_env = 2 # system E : environment system
qnum_whole = qnum_tar + qnum_ref + qnum_env
# two random states in system A+R+E (A+R:set randomly, E:set |0> initialy)
de1_whole = random_densop(qnum_tar,qnum_ref,qnum_env)
de2_whole = random_densop(qnum_tar,qnum_ref,qnum_env)
# two states in system A (trace out R+E)
de1_ini = de1_whole.partial(id=list(range(qnum_tar)))
de2_ini = de2_whole.partial(id=list(range(qnum_tar)))
# trace distance for initial states
dis_ini = de1_ini.distance(de2_ini)
# unitary transformation for whole system
U = random_unitary(qnum_whole)
de1_whole.apply(U)
de2_whole.apply(U)
# two states in system A (trace out R+E)
de1_fin = de1_whole.partial(id=list(range(qnum_tar)))
de2_fin = de2_whole.partial(id=list(range(qnum_tar)))
# trace distance for final states
dis_fin = de1_fin.distance(de2_fin)
# result
print("* trace distance(ini) =", dis_ini)
print("* trace distance(fin) =", dis_fin)
if dis_ini >= dis_fin:
print("OK!")
else:
print("NG!")
# free memory
de1_whole.free()
de2_whole.free()
de1_ini.free()
de2_ini.free()
de1_fin.free()
de2_fin.free()
前回の記事のフィデリティを計算するメソッドを、トレース距離を計算するメソッドdistanceに置き換えただけなので、詳細説明は省略します。
実行結果は以下です。
* trace distance(ini) = 0.5235122473095817
* trace distance(fin) = 0.26385081827631107
OK!
というわけで、トレース距離が収縮していることがわかります。何度も実行しましたが必ず収縮（減少）しました。
おわりに
式の展開や証明を丁寧にやったため、今回の記事は、かなり長くなってしまいました。参考文献では、要所要所でポイントになる部分の証明はなされていますが、省略されている箇所も多いです。大事と思うところは、自分で行間を埋めるようにしているのですが、今回たまたまその分量が多くなりました。ほぼ、線形代数の基本事項なのですが、とても勉強になりました（大学時代にやっていたはずなのですが、、汗）。
さて、次回の予定ですが、いまのところ「エントロピー」のお話かなーと思っています。量子情報理論の基本というからには、これは絶対に外せない話題です。というわけで、ようやく情報理論っぽい領域に突入です。
以上