Trotter の積公式で遊んでみる

目的

Trotter の積公式 の記事でごちゃごちゃ書いたことについて検証しておきたい。

復習

行列 $A$ と $B$ に対して

$$\begin{align*} \exp \left(\frac{A + B}{n}\right) \simeq \exp \left(\frac{A}{n}\right) \exp \left(\frac{B}{n}\right) \tag{1} \end{align*}$$

が成り立ち、誤差は $\mathcal{O}\left(1/n^2\right)$ で見積もれるということであった。

計算してみよう

非可換で簡単な行列として、量子力学や量子計算でよく見かける Pauli 行列を使ってみる。

$$\begin{align*} X &= \begin{pmatrix} 0 & 1 \\ 1 & 0 \end{pmatrix} \\ Z &= \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix} \\ \end{align*}$$

実際、 $XZ = \begin{pmatrix} 0 & -1 \\ 1 & 0\end{pmatrix} \neq \begin{pmatrix} 0 & 1 \\ -1 & 0\end{pmatrix} = ZX$ でこれらの行列は非可換であるが、(1) が成立するということらしい。早速なので、NumPy で試してみよう。

行列の指数関数を計算するために SciPy の expm を使うことにする。

import numpy as np
from scipy.linalg import expm

そして

X = np.array([
    [0, 1],
    [1, 0]
])
Y = np.array([
    [1, 0],
    [0, -1]
])

行列を定義する。誤差を計算するために行列のノルムを使うのだが、今回は $A = (a_{ij})$ に対して

$$\begin{align*} \| A \| = \left(\sum_{i,j} |a_{ij}|^2\right)^{1/2} \end{align*}$$

を用いる。^[1]

以下のようなコードで行列の指数関数を求め、差について $\| A - B\|^2 \times n^2$ の形で出力する。理屈通りなら、何かしら定数っぽいものが出てくることが期待できる。

for n in [3, 5, 10, 50, 100, 500]:
    A = expm((X+Y)/n)
    B = expm(X/n)@expm(Y/n)
    print(A)
    print(B)
    print('diff * n^2:', np.sqrt(np.sum((A-B)**2)) * n**2)
    print('-'*30)

実行すると以下のようになる。

[[1.45900094 0.34581691]
 [0.34581691 0.76736711]]
[[1.47386702 0.24329144]
 [0.47386702 0.75670856]]
diff * n^2: 1.4854883068493039
------------------------------
[[1.24294473 0.20267735]
 [0.20267735 0.83759003]]
[[1.24591235 0.16483998]
 [0.24591235 0.83516002]]
diff * n^2: 1.4395406249247398
------------------------------
[[1.11035034 0.10033367]
 [0.10033367 0.90968301]]
[[1.11070138 0.09063462]
 [0.11070138 0.90936538]]
diff * n^2: 1.4205105293677018
------------------------------
[[1.02040269 0.02000267]
 [0.02000267 0.98039736]]
[[1.02040539 0.01960528]
 [0.02040539 0.98039472]]
diff * n^2: 1.414464996624932
------------------------------
[[1.01010034 0.01000033]
 [0.01000033 0.99009967]]
[[1.01010067 0.00990066]
 [0.01010067 0.99009934]]
diff * n^2: 1.4142764174658962
------------------------------
[[1.002004 0.002   ]
 [0.002    0.998004]]
[[1.00200401 0.00199601]
 [0.00200401 0.99800399]]
diff * n^2: 1.4142160765324245

のような結果になる。
概ね $\| \exp \left(\frac{X + Y}{n}\right) - \exp \left(\frac{X}{n}\right) \exp \left(\frac{Y}{n}\right) \|^2 \approx \frac{1.4}{n^2}$ になっているようだ。

まとめ

単に SciPy で特定のケースについて計算してみただけだが、理論に符合する結果が出てくるのは面白い。
また実は内心本当に本当かな？と疑心暗鬼な部分があったので、何かしら確信を強められる計算結果が得られたのも良かった。

補足

$A = -iX$, $B = -iZ$ として、$\varepsilon = 1/n$ と書くと、(1) は

$$\begin{align*} \exp (-i \varepsilon X) \exp (-i \varepsilon Z) \simeq \exp ( -i \varepsilon (X+Z) ) \simeq \exp (-i \varepsilon Z) \exp (-i \varepsilon X) \end{align*}$$

という形になる。文献 [NC2] によると $R_x(\theta) = \exp (-i\theta X)$ と $R_z(\theta) = \exp (-i\theta Z)$ はそれぞれ X 軸, Z 軸周りの $\theta$ の回転ゲートということになるが、$\varepsilon$ が十分に小さい時には、この無限小の回転は $\varepsilon$ の 2 乗以降の項を無視してあたかも可換であるかのように見ることができるということを示唆していると思われる。実際、$\varepsilon$ の 1 乗の項までは一致している。

今回これ以上掘り下げるつもりもないが、同じ式を色んな言葉や状況で言い直してみるとそれだけでも味わいが増えて面白いと思う。

文献

[NC2] M. A. Nielsen, I. L. Chuang, 量子コンピュータと量子通信II, p.8, オーム社, 2005

脚注

1. 有限次元ベクトル空間上のノルムは互いに同値なので、オーダーを見る上で使いやすいやつを選んだ。