$$
\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を使って、そのプロセスの一例を再現したいと思います。
参考にさせていただいたのは、以下の記事・文献です。
一般化測定の理論
測定のプロセス
まず、何らかの対象物を測定することをイメージしてください。全体の系は、測定したい対象物を表す「注目系(S)」と、測定装置を含む「環境系(E)」とに分けられます。いま量子的な測定について考えていきたいので、「注目系(S)」も「環境系(E)」も量子系、つまり量子状態（＝密度演算子）で系の状態が完全に表現できるものとします。
このとき、測定のプロセスは、
- [STEP.1] 測定装置と対象物が独立して存在
- [STEP.2] 測定装置を対象物に当てる
- [STEP.3] 測定装置の目盛りを読み取る
- [STEP.4] 測定装置を対象物から外す
のように推移します。これを密度演算子の言葉を使って、どのように記述できるかを考えてみます。
[STEP.1]
注目系(S)の初期密度演算子を$\rho_{ini}$、環境系(E)の初期状態を$\ket{0}_E$とおくと、ただ単に、２つの系がそこに存在しているだけなので、以下の積状態で記述できます。
\rho_{ini} \otimes \ket{0}_E \bra{0}_E \tag{1}
[STEP.2]
両者間に相互作用が発生します。全体としての量子系の時間発展はユニタリ変換によって記述されるので、相互作用の結果、以下のように全体系が変化します。
\rho_{ini} \otimes \ket{0}_E \bra{0}_E \rightarrow U (\rho_{ini} \otimes \ket{0}_E \bra{0}_E) U^{\dagger} \tag{2}
このとき、注目系(S)は、どうなっているかというと、環境系(E)に関する部分トレースをとれば良いのでした。やってみます。
\begin{align}
Tr_E \{ U (\rho_{ini} \otimes \ket{0}_E \bra{0}_E) U^{\dagger} \} &= \sum_{j} \bra{j}_E U (\rho_{ini} \otimes \ket{0}_E \bra{0}_E) U^{\dagger} \ket{j}_E \\
&= \sum_{j} \bra{j}_E U \ket{0}_E \rho_{ini} \bra{0}_E U^{\dagger} \ket{j}_E \\
&= \sum_{j} M_{j} \rho_{ini} M_{j}^{\dagger} \tag{3}
\end{align}
ここで、
M_{j} \equiv \bra{j}_E U \ket{0}_E \tag{4}
と定義しました。これを「Kraus演算子」と呼び、式(3)の最終行の表現を「Kraus表現」（または、operator-sum representation）と呼びます。部分系の密度演算子の時間発展は、このような「Kraus表現」で記述することができます。
ちなみに、Kraus演算子は以下のように「完全性関係」を満たします。
\begin{align}
\sum_{j} M_{j}^{\dagger} M_{j}
&= \sum_{j} \bra{0}_E U^{\dagger} \ket{j}_E \bra{j}_E U \ket{0}_E \\
&= \bra{0}_E U^{\dagger} U \ket{0}_E \\
&= \bra{0}_E I_S \otimes I_E \ket{0}_E \\
&= I_S \tag{8}
\end{align}
また、Kraus演算子の完全性から、逆にそれを実現する全体系の演算子Uがユニタリでなければならないということも証明できます（「量子情報工学」(p.104)）。
[STEP.3]
測定装置の読み取りです。これは、環境系（例えば測定者の認識系の中）で波束の収縮が起きて、固有状態$\{ \ket{j}_E\}$の重ね合わせから１つの固有状態$\ket{n}_E$が選び取られたということに相当します。つまり、環境系(E)での射影演算です。環境系をある種のブラックボックス的なオブザーバブルと考えて、とにかく１つの固有状態に収縮するのだから、射影でいいだろ的な感じもしなくもないですが、否定する理由もないですし、ここでは、これで良いことにします。このあたり突き詰めて考えると、本当のところ、よくわからなくなってくるのですが（例えば、測定者を観測する測定者がいる場合とか）、とりあえず、今回は深みは避けて問題放置しておきます（汗）。
で、この射影演算を数式表現すると、環境系(E)で部分トレースをとる代わりに、$\bra{n}_E, \space \ket{n}_E$ではさめば良いです。やってみます。
\begin{align}
& \bra{n}_E U (\rho_{ini} \otimes \ket{0}_E \bra{0}_E) U^{\dagger} \ket{n}_E \\
&= \bra{n}_E U \ket{0}_E \rho_{ini} \bra{0}_E U^{\dagger} \ket{n}_E \\
&= M_{n} \rho_{ini} M_{n}^{\dagger} \tag{5}
\end{align}
これが注目系(S)での測定後の状態になります（ただし、正規化前）。
ここで、式(5)の形や完全性が成り立つことから、このKraus演算子というのは射影演算子にとても似ているということにご注意ください。$M_{j}$を射影演算子$P_{j}$に置き換えると、射影測定そのものの話になりますので。ということで、射影測定の射影演算子を一般化したものが、Kraus演算子である、とも言えそうです。
[STEP.4]
測定後の注目系(S)と環境系(E)の積状態になります、という一言で終了です。
というわけで、改めて各測定プロセスを、注目系(S)と環境系(E)の言葉で記載すると、以下のようになります。
- [STEP.1] 注目系(S)と環境系(E)との間に相互作用なし
- [STEP.2] 注目系(S)と環境系(E)が相互作用
- [STEP.3] 環境系(E)で射影測定
- [STEP.4] 注目系(S)と環境系(E)との間に相互作用なし
ポイントは、[STEP.2]で、部分系に注目するとその時間発展はKraus表現で表現できること。それから、[STEP.3]で、測定というのは環境系での射影演算にほかならないということです。覚えておきましょう。
確率と状態
さて、量子系の「測定」というからには、それが実現される確率はどうなの？ということに答えないといけないので、以下で答えます。
いま、初期密度演算子$\rho_{ini}$を以下のように書くことにします。
\rho_{ini} = \sum_{i} q_{i} \ket{\psi_i}_S \bra{\psi_i}_S \tag{6}
ここで、ノルム１の$\ket{\psi_i}$に対して$M_n$を演算した後のノルムが、$\ket{\psi_i}$に対する測定値nの確率です。つまり、
\bra{\psi_i}_S M_{n}^{\dagger} M_{n} \ket{\psi_i}_S = Tr(M_{n}^{\dagger} M_{n} \ket{\psi_i}_S \bra{\psi_i}_S) \tag{7}
です。初期密度演算子は、量子状態アンサンブル$\{ q_i, \ket{\psi_i}\}$のことなので、測定値nが実現される確率は、iについての和をとって、以下のようになります。
\begin{align}
p_n &= \sum_{i} q_i Tr(M_{n}^{\dagger} M_{n} \ket{\psi_i}_S \bra{\psi_i}_S) \\
&= Tr(M_{n}^{\dagger} M_{n} \rho_{ini}) \tag{8}
\end{align}
測定後の状態は、式(5)だったのですが、それは正規化前の表式でした。正しくは、式(8)の確率で割り算する必要がありますので、
\rho_{fin} = \frac{M_{n} \rho_{ini} M_{n}^{\dagger}}{Tr(M_{n}^{\dagger} M_{n} \rho_{ini})} \tag{9}
となります。
POVM測定
E_{j} \equiv M_{j}^{\dagger} M_{j} \tag{10}
のように定義される演算子を使うと、式(8)の確率は、
p_n = Tr(E_n \rho_{ini}) \tag{11}
のように書き直せます。また、Kraus演算子の完全性から、
\sum_{j} E_{j} = I \tag{12}
が成り立ちます。このような$E_j$による一般化測定の定式化もなされていて、「POVM測定(positive operator valued measure)」と呼ばれています。
一般化測定の例
参考文献の「量子情報工学（p.105）」に掲載されている例を使って、射影測定ではない、一般化測定というものをもう少し具体的にイメージしてみたいと思います。
$\ket{\phi_0}=\ket{0}$と$\ket{\phi_1}=\cos \theta \ket{0} + \sin \theta \ket{1}$という状態をもった粒子が飛んできて、どっちが来たかを判別する測定を考えます（$-\pi/2 \leq \theta \leq \pi/2$）。物理的には、$\ket{\phi_0}$は水平偏光したフォトンで$\ket{\phi_1}$は斜め直線偏光したフォトンと思っておけば良いです。
これに対して、射影演算子$P_0 = \ket{0}\bra{0},P_1=\ket{1}\bra{1}$で射影測定するのは良い選択ではありません。この測定装置でわかるのは、飛んできた粒子に$\ket{0}$成分があれば測定値が0になる可能性があり、$\ket{1}$成分があれば測定値が1になる可能性があるということです。この測定を繰り返せば、$\ket{0}$の割合と$\ket{1}$の割合はわかると思いますが、$\ket{\phi_0}$と$\ket{\phi_1}$を判別することには使えません。
そこで、POVMの登場です。やや天下りですが、以下のような演算子$E_0,E_1,E_2$を定義し、これをPOVMとして測定することを考えます。
\begin{align}
E_0 &= \frac{1}{1+\cos\theta}(\sin\theta\ket{0}-\cos\theta\ket{1})(\sin\theta\bra{0}-\cos\theta\bra{1}) \\
E_1 &= \frac{1}{1+\cos\theta} \ket{1}\bra{1} \\
E_2 &= 1-E_0-E_1 \tag{16}
\end{align}
実際にこのような測定装置を物理的にどうやって実現するかという話はちょっと置いておいて、とにかく飛んできた粒子を測定すると、0か1か2という値を出力する装置であり、測定と同時に各々の値に応じた演算子$E_0,E_1,E_2$で量子状態が変化しています、と思っておいてください。この演算子の中身を見てみると、各々の役割がわかります。$E_0$は$\ket{\phi_1}$と直交する偏光を検出します。$E_1$は$\ket{\phi_0}$と直交、すなわち垂直偏光を検出します。$E_2$はどちらでもなかったということを表しています。
したがって、測定値が0だったとすると、少なくとも$\ket{\phi_1}$ではなかった、すなわち、$\ket{\phi_0}$だったということがわかります。測定値が1だったとすると、少なくとも$\ket{\phi_0}$ではなかった、すなわち、$\ket{\phi_1}$だったということがわかります。測定値が2だったとすると、これはどちらとも判定できなかったということなので、判別不能と考えます。
というわけで、判別不能という場合もあるのですが、間違って判別することはありません。このように、一般化測定を考えると、非直交状態でも誤りのない判別（USD:unambiguous state discremination）が可能になります。ちなみに、$\theta=\pi/2$の場合、容易にわかるように、$\ket{0},\ket{1}$を判別する単なる射影測定になります（$E_0,E_1$は射影演算子、$E_2$は0なので意味がなくなります）。
一般化測定のシミュレーション
それでは、たったいま説明したPOVM測定の例を、量子計算シミュレータqlazyで再現してみたいと思います。どっちかの粒子が飛んできたときに、どのくらいの確率で判別が成功するのかを計算してみます。これをやるため、qlazyに一般化測定する機能を追加しました（v0.0.24）。
全体のPythonコードを示します。
【2021.9.5追記】qlazy最新版でのソースコードはここに置いてあります。
import math
import numpy as np
from qlazypy import QState, DensOp
def make_povm(theta=0.0):
s = math.sin(theta * math.pi)
c = math.cos(theta * math.pi)
f = 1.0/(1.0 + c)
E0 = f * np.array([[s*s, -c*s], [-c*s, c*c]])
E1 = f * np.array([[0, 0], [0, 1]])
E2 = np.eye(2) - E0 - E1
return (E0, E1, E2)
if __name__ == '__main__':
theta = 0.5 # -0.5 <= theta <= 0.5
qs = [QState(1), QState(1).ry(0,phase=theta*2.0)]
de = [DensOp(qstate=[qs[0]], prob=[1.0]), DensOp(qstate=[qs[1]], prob=[1.0])]
povm = make_povm(theta)
print("theta = {0:}*PI".format(theta))
print(" [0] [1] [?]")
for i in range(2):
prob = de[i].probability(povm=povm)
print("[{0:}] {1:.3f} {2:.3f} {3:.3f}".format(i, prob[0],prob[1],prob[2]))
de[i].free()
qs[i].free()
何をやっているか、説明します。
theta = 0.5 # -0.5 <= theta <= 0.5
$\theta$の値を設定します。角度の単位は$\pi$ラジアンです。なので、theta=0.5は$0.5\pi$ラジアンを表します。
qs = [QState(1), QState(1).ry(0,phase=theta*2.0)]
量子状態を表すQStateクラスで、２種類の状態を作成してリスト化しています。QState(1)で１つの粒子が$\ket{0}$となっている状態を作成します。QState(1).ry(0,phase=theta*2.0)で１つの粒子が$\ket{0}$となっている状態を作成し、それにY軸周りの回転ゲートで$\ket{\phi_1}$を作成するようしています。ということで、qs[0]は$\ket{\phi_0}$、qs[1]は$\ket{\phi_1}$を表す形になります。
de = [DensOp(qstate=[qs[0]], prob=[1.0]), DensOp(qstate=[qs[1]], prob=[1.0])]
DensOpクラスで、その量子状態各々に対する密度演算子を作成します。各々１個の純粋状態に対する密度演算子なので、指定する量子状態は１つ、確率も1.0と指定します。これで、$\ket{\phi_0}$に対する密度演算子がde[0]、$\ket{\phi_1}$に対する密度演算子がde[1]と定義されることになります。
povm = make_povm(theta)
POVM測定の演算子をnumpy行列のタプルの形で変数povmに格納します。詳細は関数定義の部分を見てください。３つの演算子$E_0,E_1,E_2$を作成しています。
print("theta = {0:}*PI".format(theta))
print(" [0] [1] [?]")
for i in range(2):
prob = de[i].probability(povm=povm)
print("[{0:}] {1:.3f} {2:.3f} {3:.3f}".format(i, prob[0],prob[1],prob[2]))
print文は結果表示をわかりやすくするためのものなので、一旦無視していただき、forループの中身に注目してください。
prob = de[i].probability(povm=povm)
密度演算子に対してprobabilityメソッドで一般化測定をシミュレーションします。引数として与えるのは、先程定義したPOVM演算子に対応した、numpy行列のタプルです。返り値は各演算子に対応した確率値のリストです。いま３つの要素からなるPOVMを実行しているので、確率として３つの値（$E_0,E_1,E_2$に対応した）がリストの要素になります。
これをforループで回します。はじめに$\ket{\phi_0}$が飛んできたときに、どんな確率で判別されるかが計算され、次に$\ket{\phi_1}$が飛んできたときに、どんな確率で判別されるかが計算されて、終了です。
結果はどうなるかというと、、、
まず、$\theta=0.5\pi$の場合です。
theta = 0.5*PI
[0] [1] [?]
[0] 1.000 0.000 0.000
[1] 0.000 1.000 0.000
ここで、縦に並んでいる[0],[1]は、各々$\ket{\phi_0}$が飛んできた場合、$\ket{\phi_1}$が飛んできた場合を表しています。横に並んでいる[0],[1],[?]は、各々$\ket{\phi_0}$と判別される確率、$\ket{\phi_1}$と判別される確率、判別不能である確率を表しています。数字はその確率値です。
$\theta=0.5\pi$は両者直交している場合なので、完璧に判別できることがわかります。
次に、$\theta=0.4\pi$と$\theta=0.3\pi$の場合です。
theta = 0.4*PI
[0] [1] [?]
[0] 0.691 0.000 0.309
[1] 0.000 0.691 0.309
theta = 0.3*PI
[0] [1] [?]
[0] 0.412 0.000 0.588
[1] 0.000 0.412 0.588
ということで、角度が小さくなる（両者平行に近くなる）に従い、判別できる確率が小さくなることがわかります。が、決して判別ミスすることはないです。つまり、$\ket{\phi_0}$が来たのに$\ket{\phi_1}$であると判別することはありません、ということもわかります。
最後に、$\theta=0.0$の場合です。
theta = 0.0*PI
[0] [1] [?]
[0] 0.000 0.000 1.000
[1] 0.000 0.000 1.000
完璧に判別不能となりました。両者平行なので当然の結果です。
おわりに
注目している対象物の現実的な測定は、実はその周りの環境系で射影測定した結果である、というお話でした。一言で言うととてもシンプルなのですが、この理解を進める中で、部分系の時間発展は完全性を満たす「Kraus演算子」を使った「Kraus表現」で記述されるという知見も得ました。一方で、密度演算子は必ず半正定値でトレースが１になるのでした。数学的には両者微妙に違うことを言っているような気がするのですが、実は同じことのようです。「CPTPマップ」というものを勉強すれば、理解できそうなので、次はこれかなー、と思っています。
以上