Lindblad 方程式の導出

---

1. Lindblad 方程式

Lindblad 方程式とは、デコヒーレンスのある系で密度演算子の時間発展を記述する方程式で、以下のような形である。 \[\frac{d\rho}{dt} = -i[H,\rho(t)] + \sum_\mu \sum_\nu c_{\mu\nu}\left(F_\mu\rho F_\nu^\dagger - \frac{1}{2}\left\{F_\mu^\dagger F_\nu , \rho\right\}\right) \] 今回は、色々な近似を使いながら、シュレディンガー方程式からこの Lindblad 方程式の導出をしてみようと思う。

---

---

2. 考え方

デコヒーレンスを理論的に考えるときには、注目している量子系が、注目の外にある「環境」と相互作用することによって、量子状態\(\ket{\psi}\)が壊れてしまうという考え方をする。

考えている系を \(S\) (system) と書き、環境は \(B\) (bath) と書く。\(S,B\)をあわせた系のハミルトニアンは一般に、 \[H = H_S + H_B + H_{SB} \tag{1}\] という形をしているだろう。\(H_S, H_B\)はそれぞれ注目している系と環境のハミルトニアン、\(H_{SB}\)が系と環境の相互作用を表すハミルトニアンとする。

普通の量子力学では、暗に\(H_B, H_{SB}\)が完全に無視できる、という仮定を入れている。そういう単純なモデル化に成功したからこそ、シュレディンガー方程式 \[\frac{d}{dt}\ket{\psi(t)}_S = H_S \ket{\psi(t)}_S\] は成功を収めたわけである。

しかし本来は、環境と系を含めた(1)のハミルトニアンにしたがって時間発展すると考えられる。つまり全体系の時間発展は \[\frac{d}{dt}\ket{\psi(t)} = H \ket{\psi(t)}\] となるはずだ。

そんな複雑なハミルトニアンで時間発展している中で、僕たちが測定できるのは注目している系 \(S\) だけである。そのような場合には、環境 \(B\) の影響がある種平均化されたような形で、\(S\)の測定値に現れる。現在のところ、これがデコヒーレンスの原因であると考えられている。

---

3. 使う近似・仮定

僕が勉強した時には、次々に色々な近似や仮定が現れてわけがわからなくなったので、先にどんな近似や仮定を使うかまとめておく。\(S\)の密度行列を\(\rho_S\)、\(B\)の密度行列を\(\rho_B\)と書き、全体系の密度行列は\(\rho\)と書くことにする。

---

4. 相互作用表示

系\(S\)と環境\(B\)全体を動かすハミルトニアンは \[H = H_S + H_B + H_{SB} \tag{1}\] であり、全体系の密度演算子は \[\frac{d\rho}{dt} = -i[H,\rho]\] という微分方程式に従う。このページでは\(\hbar=1\)とする。相互作用ハミルトニアン\(H_{SB}\)は、 \[H_{SB} = \sum_j S_j \otimes B_j\] という形をしているとしよう。\(S_j, B_j\)はそれぞれ\(S\)と\(B\)に作用するエルミート演算子とする。

このように相互作用する系を解析するには、相互作用表示 が便利だ。今回の場合、相互作用表示の密度演算子\(\rho_I(t)\)を、 \[\rho_I(t) = e^{i(H_S+H_B)t} \rho(t) e^{-i(H_S+H_B)t}\] と置くことにより、(1)のハミルトニアンのうち、相互作用項だけを残した座標系に移ることができる。上の微分方程式に代入してみればすぐに分かるように、\(\rho_I(t)\)は以下の式に従う。 \[\frac{d\rho_I}{dt} = -i [ e^{i(H_S+H_B)t}H_{SB} e^{-i(H_S+H_B)t}, \rho_I(t)] \] 式を簡単にするために、相互作用表示における相互作用ハミルトニアンを \[H_I(t) = e^{i(H_S+H_B)t}H_{SB} e^{-i(H_S+H_B)t} \tag{2}\] と書くことにする。そうすると、相互作用表示における時間発展の方程式として、 \[\frac{d\rho_I}{dt} = -i [ H_I(t), \rho_I(t)] \tag{3}\] を得る。この式をもとにLindblad方程式を導出していこう。

---

5. 導出

目標は注目している系\(S\)の密度演算子\(\rho_S\)に対する微分方程式を作ることである。\(\Tr_B (\cdot) \)を環境 \(B\) についての部分トレースとして、\(\rho_S = \Tr_B (\rho)\)と表せることに注意しておこう。

まず ダイソン級数 を導出したときと同じような式変形をする。(3)式を積分表示すると、 \[\rho_I(t) = \rho_I(0) -i \int_0^t [H_I(\tau), \rho_I(\tau)] d\tau\] となり、これを(3)式の右辺に代入してやると、次のような表式が得られる。 \[\frac{d\rho_I}{dt} = -i [ H_I(t), \rho_I(0)] - \int_0^t [H_I(t),[H_I(\tau),\rho_I(\tau)]] d\tau\] さて、さらに両辺(の\(\Tr_B (\cdot) \)をとって\(\rho_S\)に関する微分方程式に直すと、相互作用表示における系\(S\)の密度演算子を\(\rho_{SI}(t) = e^{iH_St} \rho_S(t) e^{-iH_St}\)として \[\frac{d\rho_{SI}}{dt} = -i \Tr_B \left([ H_I(t), \rho_I(0)]\right) - \int_0^t \Tr_B\left([H_I(t),[H_I(\tau),\rho_I(\tau)]]\right) d\tau\] を得る。

ここで、実は常に \[\Tr_B \left([ H_I(t), \rho_I(0)]\right) = 0\] を満たすように相互作用表示を取ることができる。これを満たすようにする具体的なやり方をここで書くと、導出の本筋からそれてしまうので、このページ参照。ともかく、一般性を失わずに、\(\rho_{SI}(t)\)は \[\frac{d\rho_{SI}}{dt} = - \int_0^t \Tr_B\left([H_I(t),[H_I(\tau),\rho_I(\tau)]]\right) d\tau\] という形の方程式を満たすといえる。

さて、ここからは次々に近似を行う。まずボルン近似によって、\(\rho_I(t) \approx \rho_{SI}(t)\otimes\rho_B(0)\)として、 \[\frac{d\rho_{SI}}{dt} = - \int_0^t \Tr_B\left([H_I(t),[H_I(\tau),\rho_{SI}(\tau)\otimes\rho_B(0)]]\right) d\tau\] さらにマルコフ近似により\(\rho_{SI}(\tau) \to \rho_{SI}(t)\)とすると \[\frac{d\rho_{SI}}{dt} = - \int_0^t \Tr_B\left([H_I(t),[H_I(\tau),\rho_{SI}(t)\otimes\rho_B(0)]]\right) d\tau \tag{4}\] を得る。これで\(\rho_{SI}(t)\)が現在の状態のみに依存して変化する方程式を得ることができた。

ボルン近似はともかく、個人的にはマルコフ近似は「近似」といっていいものなのかわからない。特に理由もなさそうなので、どちらかというと「仮定」に近いものなのではないだろうか。まあ過去のすべての状態に依存するような方程式は、解析するのが面倒なのは確かなのでそうしてしまうのだろう。そんなことよりもむしろ、ページの一番最初にあげたLindblad方程式は、\(\frac{d\rho_S}{dt} = \mathcal{L}\rho_S(t)\) (\(\mathcal{L}\)は線形演算子) という形の方程式の中で、密度行列のトレースと半正定値性を保存するようなものの一般系を与えていることが知られているので、それに合わせこむように適当に理由づけているのかもしれない。

次に\(H_{SB} = \sum_j S_j \otimes B_j\)を代入してもう少し具体的な形を考えてみよう。\(H_I(t)\)は\(H_{SB}\)を相互作用表示したものだったので、相互作用表示における\(S_j, B_j\)を \(S_j(t) = e^{iH_St}S_je^{-iH_St}\)、\(B_j(t) = e^{iH_St}B_je^{-iH_St}\)と書くことにすると、 \[H_I(t) = \sum_j S_j(t) \otimes B_j(t) \tag{5}\] である。これを(4)に代入して、計算を進めていくと \begin{align} \frac{d\rho_{SI}}{dt} &= - \int_0^t \Tr_B\left(\left[\sum_j S_j(t) \otimes B_j(t),\left[\sum_k S_k(\tau) \otimes B_k(\tau),\rho_{SI}(t)\otimes\rho_B(0)\right]\right]\right) d\tau \\ &= - \sum_j \sum_k\int_0^t \Tr_B\left(\left[S_j(t) \otimes B_j(t),\left[S_k(\tau) \otimes B_k(\tau),\rho_{SI}(t)\otimes\rho_B(0)\right]\right]\right) d\tau \\ &= - \sum_j \sum_k\int_0^t \Tr_B\left(\left[S_j(t) \otimes B_j(t), S_k(\tau)\rho_{SI}(t) \otimes B_k(\tau)\rho_B(0) - \rho_{SI}(t)S_k(\tau) \otimes \rho_B(0)B_k(\tau)\right]\right) d\tau \\\\ &= - \sum_j \sum_k\int_0^t \Tr_B\left(S_j(t)S_k(\tau)\rho_{SI}(t) \otimes B_j(t)B_k(\tau)\rho_B(0) - S_k(\tau)\rho_{SI}(t)S_j(t) \otimes B_k(\tau)\rho_B(0)B_j(t)\right. \\ &~~~ \left.- S_j(t)\rho_{SI}(t)S_k(\tau) \otimes B_j(t)\rho_B(0)B_k(\tau) + \rho_{SI}(t)S_k(\tau)S_j(t) \otimes \rho_B(0)B_k(\tau)B_j(t)\right) d\tau \\\\ &= - \sum_j \sum_k\int_0^t \left(S_j(t)S_k(\tau)\rho_{SI}(t) \Tr(B_j(t)B_k(\tau)\rho_B(0)) - S_k(\tau)\rho_{SI}(t)S_j(t) \Tr (B_k(\tau)\rho_B(0)B_j(t))\right. \\ &~~~ \left.- S_j(t)\rho_{SI}(t)S_k(\tau) \Tr (B_j(t)\rho_B(0)B_k(\tau)) + \rho_{SI}(t)S_k(\tau)S_j(t) \Tr (\rho_B(0)B_k(\tau)B_j(t))\right) d\tau \\ \end{align} だんだん式が複雑になってきたので、ここらで記号を導入しておこう。 \[G_{jk} (t,\tau) := \Tr(B_j(t)B_k(\tau)\rho_B(0)) \tag{6}\] と書くことにする。この関数 \(G_{jk} (t,\tau)\) は

相関関数

と呼ばれる。 相関関数とは、その名のとおり\(B_j(t)\)と\(B_k(t)\)の「相関」を測る関数である。今回の場合に落とし込んでもう少し具体的にすると、時刻\((t,\tau)\)における相関関数 \(G_{jk}(t,\tau)\) が正ならば、\(B_j(t)\)と\(B_k(\tau)\)は同じ符号を取ることが「多い」と考えられ、負ならば逆符号を取ることが多いと考えられる。逆に相関関数が 0 であるということは、その2つには全く相関が無く、片方からもう一方の様子を窺い知ることができないような状況を表す。直感的に、\(G_{jk}(t,\tau)\) はその系が十分複雑で、時刻差 \(|t-\tau|\) が十分大きいときにほとんど\(G_{jk}(t,\tau)\approx 0\)となることが予想されるだろう。複雑な系なら、\(B_j(t)\)と\(B_k(t)\)が長い時間関係を保っているとは考えにくい。

---

上の式をみると、(6)で定義した\(G_{jk} (t,\tau)\)だけではすべての項を表せないように見えるが、実はできる。トレースの性質と\(\{B_j\}\)のエルミート性を使うと、 \[\Tr (B_k(\tau)\rho_B(0)B_j(t)) = \Tr (B_j(t)B_k(\tau)\rho_B(0)) = G_{jk} (t,\tau)\] や \begin{align} \Tr (\rho_B(0)B_k(\tau)B_j(t)) &= \left(\Tr \left((\rho_B(0)B_k(\tau)B_j(t))^\dagger\right)\right)^* \\ &= \left(\Tr \left(B_j^\dagger(t)B_k^\dagger(\tau)\rho_B^\dagger(0)\right)\right)^*\\ &= \left(\Tr \left(B_j(t)B_k(\tau)\rho_B(0)\right)\right)^*\\ &= G_{jk}^* (t,\tau) \end{align} が得られる。したがって \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k\int_0^t d\tau\left[G_{jk}(t,\tau)\left(S_j(t)S_k(\tau)\rho_{SI}(t) - S_k(\tau)\rho_{SI}(t)S_j(t) \right) + H.c.\right]\tag{7} \end{align} と書けることがわかるだろう。この方程式は

Redfield 方程式

と呼ばれることもあるようだ。この式において\(\tau\)は\(t\)以前の時刻からの寄与をすべて集める役割をしている。そのことをさらに明確にするために、「どれだけ過去の時刻にさかのぼるか」を表すパラメータ\(s\)を\(\tau = t-s\)によって導入しよう。そして積分を\(s\)に関するものに書き換えると \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k\int_0^t ds\left[G_{jk}(t,t-s)\left(S_j(t)S_k(t-s)\rho_{SI}(t) - S_k(t-s)\rho_{SI}(t)S_j(t) \right) + H.c.\right] \end{align} を得る。このように変形してやることで、環境の相関時間が十分に小さいという近似が使いやすくなっている。実験的に観測する時間\(t\)が、相関関数が減衰する時間\(\tau_B\)よりも十分大きいならば、積分範囲を\(0\to t\)から\(0 \to \infty\)としてもほとんど誤差はないはずである。そこで \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k\int_0^\infty ds\left[G_{jk}(t,t-s)\left(S_j(t)S_k(t-s)\rho_{SI}(t) - S_k(t-s)\rho_{SI}(t)S_j(t) \right) + H.c.\right] \end{align} としてしまおう。さらに、相関関数が時間変化しないという仮定を使って \[G_{jk}(s) := G_{jk}(s,0)\] と書くことにし、 \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k\int_0^\infty ds\left[G_{jk}(s)\left(S_j(t)S_k(t-s)\rho_{SI}(t) - S_k(t-s)\rho_{SI}(t)S_j(t) \right) + H.c.\right]\tag{8} \end{align} を得る。

H.c.はHermitian conjugateのことで、その前の項のエルミート共役を表す。

---

最後に Secular Approximationを施してやればLindblad方程式が出てくる。Secular Approximationとは、高速に振動する項は落としてしまう近似であったので、まずは\(\{S_j(t)\}\)を振動数ごとに分解してやらないといけない。

\(S_j(t)\)は \[S_j(t) = e^{iH_St}S_je^{-iH_St}\] と定義されていた。振動数ごとに分解するために、\(H_S\)の固有状態\(\ket{\omega}\) (それぞれの固有値は\(\omega\)とする) によって\(S_j\)を行列表記してやろう。 \[S_j = \sum_{\omega,\omega'}\ket{\omega}\bra{\omega}S_j\ket{\omega'}\bra{\omega'}\] とすれば \begin{align} S_j(t) &= \sum_{\omega,\omega'} e^{i\omega t}\ket{\omega}\bra{\omega}S_j\ket{\omega'}\bra{\omega'}e^{-i\omega' t} \\ &= \sum_{\omega,\omega'} e^{-i(\omega' - \omega) t}\ket{\omega}\bra{\omega}S_j\ket{\omega'}\bra{\omega'} \end{align} ここから、足し算の順番を少し変えて、ある\(\Omega\) について \(\omega' - \omega = \Omega\)となる部分のみを先に足しあげて、あとで\(\Omega\)に関する和をとれば、\(S_j(t)\)を各周波数\(\Omega\)の成分に分解できる。すなわち \begin{align} S_j(t) &= \sum_{\Omega}\sum_{\omega' - \omega = \Omega} e^{-i(\omega' - \omega) t}\ket{\omega}\bra{\omega}S_j\ket{\omega'}\bra{\omega'}\\ &= \sum_{\Omega}e^{-i\Omega t}\sum_{\omega' - \omega = \Omega} \ket{\omega}\bra{\omega}S_j\ket{\omega'}\bra{\omega'} \tag{9} \end{align} を得て、これは確かに、2つ目の和の部分 \begin{align} S_j(\Omega) &=\sum_{\omega' - \omega = \Omega} \ket{\omega}\bra{\omega}S_j\ket{\omega'}\bra{\omega'} \end{align} を係数として周波数分解された形になっている。ちなみに、\(S_j^\dagger(\Omega) = S_j(-\Omega)\)に注意しよう。

さて、この(9)式を(8)へ代入して \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k \sum_\Omega \sum_{\Omega'}\int_0^\infty ds\left[G_{jk}(s)\left( e^{-i(\Omega+\Omega')t} e^{i\Omega s}S_j(\Omega')S_k(\Omega)\rho_{SI}(t) - e^{-i(\Omega+\Omega')t} e^{i\Omega s} S_k(\Omega)\rho_{SI}(t)S_j(\Omega') \right) + H.c.\right]\\ &= - \sum_j \sum_k \sum_\Omega \sum_{\Omega'}\int_0^\infty ds\left[e^{-i(\Omega+\Omega')t} e^{i\Omega s}G_{jk}(s)\left( S_j(\Omega')S_k(\Omega)\rho_{SI}(t) -S_k(\Omega)\rho_{SI}(t)S_j(\Omega') \right) + H.c.\right] \end{align} \(s\)の積分を \[\Gamma_{jk}(\Omega) = \int_0^\infty ds e^{i\Omega s}G_{jk}(s)\] と置くと、(この相関関数をフーリエ変換した関数は

スペクトル密度

とも呼ばれる) \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k \sum_\Omega \sum_{\Omega'}\Gamma_{jk}(\Omega)\left[e^{-i(\Omega+\Omega')t} \left( S_j(\Omega')S_k(\Omega)\rho_{SI}(t) -S_k(\Omega)\rho_{SI}(t)S_j(\Omega') \right) + H.c.\right] \end{align} を得る。最後にSecular Approximationとして\(\Omega+\Omega' = 0\)となる部分だけを取り出すと \begin{align} \frac{d\rho_{SI}}{dt} &= - \sum_j \sum_k \sum_\Omega \Gamma_{jk}(\Omega)\left[ \left( S_j(-\Omega)S_k(\Omega)\rho_{SI}(t) -S_k(\Omega)\rho_{SI}(t)S_j(-\Omega) \right) + H.c.\right] \\ &= - \sum_j \sum_k \sum_\Omega \Gamma_{jk}(\Omega)\left[ \left( S_j^\dagger(\Omega)S_k(\Omega)\rho_{SI}(t) -S_k(\Omega)\rho_{SI}(t)S_j^\dagger(\Omega) \right) + H.c.\right] \\ \end{align}

---

かなり近づいた。もう一歩だ。\(\Gamma_{jk}(\omega)\)を実部と虚部に分解し \[\Gamma_{jk}(\omega) = \frac{1}{2} \left(J_{jk}(\omega) + iI_{jk}(\omega)\right)\] として代入し、丁寧に計算すると、Lindblad方程式が得られる。 \begin{align} \frac{d\rho_{SI}}{dt} &= \sum_\omega \sum_j \sum_k \left[-\frac{i}{2}[I_{jk}(\omega)S_j^\dagger(\omega)S_k(\omega), \rho_{SI}(t)] + I_{jk}(\omega)\left(S_k(\omega)\rho_{SI}(t)S_j^\dagger(\omega) - \frac{1}{2}\left\{S_j^\dagger(\omega)S_k(\omega),\rho_{SI}(t)\right\} \right)\right] \\ \end{align} 実はここで現れたハミルトニアン的な部分 \[H_{LS} := \sum_\omega \sum_j \sum_k \frac{1}{2}I_{jk}(\omega)S_j^\dagger(\omega)S_k(\omega) \] は

ラムシフト

と呼ばれるものである。そこでこれを用いて書き直すと \begin{align} \frac{d\rho_{SI}}{dt} &= -i[H_{LS}, \rho_{SI}(t)] + \sum_\omega \sum_j \sum_k \left[I_{jk}(\omega)\left(S_k(\omega)\rho_{SI}(t)S_j^\dagger(\omega) - \frac{1}{2}\left\{S_j^\dagger(\omega)S_k(\omega),\rho_{SI}(t)\right\} \right)\right] \\ \end{align} と少しかっこよくまとまる。これは相互作用表示だが、もちろん元の座標系に戻ることもできるだろう。演習問題にちょうどいいかもしれない。

にしても近似使いすぎだよなあ。研究者たちが納得しないのも納得。こういう開放量子系の話は今ホットトピックのようだ。