$$
\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}}
$$
前置き
いつもとタイトルの雰囲気を変えてますが、今回も光量子計算の話題です。
光の量子状態に演算子(量子ゲート)を作用させる計算について、こちらの記事ではガウス状態間の変換に限りガウス分布の平均、分散に対する比較的単純な計算で書けることを紹介しました。
今回はガウス状態に限定されない、より一般的な形式の実装を試みます。
演算子の記述
一般的に、光量子における演算子は生成・消滅演算子を使って記述します。
今回は簡単な例としてDisplacement operatorを扱います。
$D(\alpha) = \exp(\alpha \hat{a}^{\dagger} - \alpha^* \hat{a})$
生成/消滅演算子$\hat{a}^{\dagger}, \hat{a}$
は、それぞれ光子数状態に対し以下のように作用します。
ざっくりいうと、生成演算子は光子を増やし消滅演算子は光子を減らします。
\hat{a} \ket{n} = \sqrt{n}\ket{n-1} \\
\hat{a} \ket{0} = 0 \\
\hat{a}^{\dagger} \ket{n} = \sqrt{n + 1}\ket{n + 1}
ならば量子状態は光子数基底で表すのが良さそうですね。
ここからは、光子数基底で表された量子状態に対してDisplacement operatorを作用させる計算を実装しながら追っていきます。
本来光子数基底は無限次元ですが、実装上それは不可能なので有限の光子数(=cutoff)以下で近似します。
import
これだけ使います。
import numpy as np
from scipy import special
消滅演算子
消滅演算子をk回、n光子数状態に作用させる関数を実装します。
\hat{a}^k \ket{n}
def down(order, n):
if order > n:
return [0, 0]
n_ = n - order
coeff = np.prod(np.sqrt(np.arange(n, n - order, -1)))
return [coeff, n_]
生成演算子
消滅演算子と同様、生成演算子をk回、n光子数状態に作用させる関数を実装します。
def up(order, n, cutoff = 10):
n_ = n + order
if n_ > cutoff:
return [0, 0]
coeff = np.prod(np.sqrt(np.arange(n + 1, n + 1 + order)))
return [coeff, n_]
指数関数の肩に演算子
指数関数の定義に従いべき級数に展開するのが良いです。
\exp(\alpha \hat{a}) = \sum_{n} \frac{1}{n!}(\alpha \hat{a})^n
先に実装したup(), down()を使って$\exp(\alpha \hat{a})(\sum_n C_n \ket{n})$, $\exp(\alpha \hat{a}^{\dagger})(\sum_n C_n \ket{n})$
をそれぞれ実装します。
入力配列の与え方ですが、入力光子数状態$\sum_n C_n \ket{n}$として、$\ket{i}$の係数$C_i$を入力配列(fockState)のi番目の要素に格納しています。
出力も同様です。
def exp_annihiration(fockState, alpha, cutoff):
state = np.zeros(cutoff + 1) + 0j
for j in range(len(fockState)):
for i in range(cutoff + 1):
tmp = down(i, j)
state[np.int(tmp[1])] += tmp[0] * fockState[j] / special.factorial([i]) * alpha ** i
return state
def exp_creation(fockState, alpha, cutoff):
state = np.zeros(cutoff + 1) + 0j
for j in range(len(fockState)):
for i in range(cutoff + 1):
tmp = up(i, j)
state[np.int(tmp[1])] += tmp[0] * fockState[j] / special.factorial([i]) * alpha ** i
return state
ここでDisplacement operatorをべき級数展開しようとすると気づくのですが、$D(\alpha) = \exp(\alpha \hat{a}^{\dagger} - \alpha^* \hat{a})$
のように指数関数の肩に2項以上あると項の数が増えます。
$\exp(\alpha \hat{a}^{\dagger}) \exp(-\alpha^* \hat{a})$
と変形すれば済むと考えたくなりますが、悲しいことにできません。
生成/消滅演算子の間にボソン由来の交換関係$[\hat{a}, \hat{a}^{\dagger}] = 1$が成り立つため、簡単にいうと非可換なためです。
そこで今回はoperator theorem
$e^{A+B} = e^{A}e^{B}e^{-[A, B]/2}$
を使って$D(\alpha)$を次のように"ordering"します。
D(\alpha) = \exp(-|\alpha|^2 / 2) \exp(\alpha \hat{a}^{\dagger}) \exp(-\alpha^* \hat{a})$
これでかなり実装しやすくなりました。
先に実装した関数を使うと、$D(\alpha)(\sum_n C_n \ket{n})$は以下のように実装できます。
def displacement(fockState, alpha, cutoff):
state_ = exp_annihiration(fockState, -np.conj(alpha), cutoff)
state = exp_creation(state_, alpha, cutoff)
return state * np.exp(-np.abs(alpha)**2 / 2)
簡単な例として、$D(\alpha)\ket{0} = \ket{\alpha}$はこの前明倫館で2000円で拾った本に展開済の式がありました。
$D(\alpha)\ket{0} = \exp{(-|\alpha|^{2} / 2)} \sum_n \frac{\alpha^n}{\sqrt{n!}} \ket{n}$
こちらも比較確認用に名前を変えて実装します。
def coherentState(alpha, cutoff = 10):
dim = cutoff + 1
n = np.arange(dim)
state = np.exp(- 0.5 * np.abs(alpha)**2) / np.sqrt(special.factorial([n])) * alpha**n
return state
さあ、結果を比較しましょう。
displacement()には光子数状態$\ket{0}$を入力します。
cutoff = 10
fockState = np.zeros(10)
fockState[0] = 1
alpha = 1 + 1j
print(coherentState(alpha, cutoff))
print(displacement(fockState, alpha, cutoff))
結果は以下の通りです。しっかり一致してますね。
テストとしてはまだかなりザルですが。
[[ 0.36787944+0.j 0.36787944+0.36787944j 0. +0.5202601j
-0.30037231+0.30037231j -0.30037231+0.j -0.13433058-0.13433058j
0. -0.10968046j 0.04145532-0.04145532j 0.02931334+0.j
0.00977111+0.00977111j 0. +0.00617979j]]
[ 0.36787944+0.j 0.36787944+0.36787944j 0. +0.5202601j
-0.30037231+0.30037231j -0.30037231+0.j -0.13433058-0.13433058j
0. -0.10968046j 0.04145532-0.04145532j 0.02931334+0.j
0.00977111+0.00977111j 0. +0.00617979j]
原則的かつ楽観的には、他の量子ゲートも同様に計算を拡張していけば実装可能だと思います。
もっと大きく大胆に言うと、生成消滅演算子で書ける任意のハミルトニアンによる時間発展が計算できるはずです。
これで！あなたも！ボソニック量子系エンジニアになれる！多分！
この実装は上手く拡張できたら下のライブラリに実装していきます。
どうぞよろしくということで。
https://github.com/Blueqat/Photonqat
以上、読んで頂きありがとうございました。