シュレーディンガー方程式の時間発展を量子回路で表現

$$
\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}}
$$

## はじめに

ハミルトニアンが時間に依存していない場合、波動関数が時間に伴ってどう変化するか、解析的に解けたり解けなかったりするわけですが、電子スピンのような２準位の多体系については、量子回路を使って近似的に表現できるようなので、自分の脳内整理の意味も込めて、まとめてみたよ。

以下を参考にさせていただきました。

* [Quantum Native Dojo: 4-2. トロッター分解を用いた量子シミュレーション](https://github.com/qulacs/quantum-native-dojo)
* [【量子化学】制限付きVQE法(1/2)](https://blog.mdrft.com/post/1620)
* [中山茂、「量子アルゴリズム」、技報堂出版(2014)](http://gihodobooks.jp/book/3343-0.html)


## 問題設定

まず、シュレディンガー方程式です。

```math
i \frac{\partial}{\partial t} \ket{\phi(t)} = H \ket{\phi(t)}
```

ここで、$H$はハミルトニアン演算子を表します。これを解くと、

```math
\begin{align}
\ket{\phi(t)} &= exp(-i H t) \ket{\phi(0)} \\
&= U(t) \ket{\phi(0)}
\end{align}
```

となります。$U(t)=exp(-i H t)$のことを時間発展演算子と呼びます。いま、$H$として、以下のようなものを考えます。すなわち、スピンを持った$N$個の粒子同士が相互作用している系です。

```math
\begin{align}
&H = \sum_{m=0}^{M-1} P_{m} \\
&P_{m} = \sigma^{p_{m}^{0}}_{0} \otimes \sigma^{p_{m}^{1}}_{1} \cdots \otimes \sigma^{p_{m}^{N-1}}_{N-1}
\end{align}
```

ここで、$\sigma^{p}_{n} (p=0,1,2,3)$は、$n$番目の粒子に関するパウリのスピン演算子を表しており、

```math
\begin{align}
\sigma_{n}^{0} &= I_{n} =
	       \begin{pmatrix}
	       1 & 0 \\
	       0 & 1
	       \end{pmatrix} \\
\sigma_{n}^{1} &= X_{n} =
	       \begin{pmatrix}
	       0 & 1 \\
	       1 & 0
	       \end{pmatrix} \\
\sigma_{n}^{2} &= Y_{n} =
	       \begin{pmatrix}
	       0 & -i \\
	       i & 0
	       \end{pmatrix} \\
\sigma_{n}^{3} &= Z_{n} =
	       \begin{pmatrix}
	       1 & 0 \\
	       0 & -1
	       \end{pmatrix}  
\end{align}

```

で定義されるものとします。

いま一度、上の時間発展演算子を見てください。$e$の肩にかかっている$H$はエルミートなので、$exp(-i H t)$はユニタリです。ユニタリであれば、なんとか量子回路で表現できるはずです。ですが、具体的に回路に落としこもうとすると、それほど自明ではありません。が、上に示したスピン相互作用のハミルトニアンであれば、近似的に量子回路で表現できるみたいです。どうやれば、そんなことができるのかを、以下で見ていきたいと思います。

が、その前に、いくつかの定義と知っておくと便利な公式を示しておきます。


## 定義と公式

### ユニタリ行列に対するオイラーの公式

$U$をユニタリ行列とすると、

```math
e^{i\theta U} = \cos(\theta) I + i \sin(\theta) U
```

が成り立ちます(ここで、$I$は単位行列です)。ちなみに、行列版の指数関数は、普通の指数関数のテイラー展開と同じように、以下の形で定義されます。

```math
e^{A} = \sum_{k=0}^{\infty} \frac{1}{k!} A^{k}
```


### パウリのスピン行列

先程、定義を示しました。以下の式が成り立ちます。

```math
\begin{align}
&X^{\dagger} = X, \; Y^{\dagger} = Y, \; Z^{\dagger} = Z \\
&X^{2} = I, \; Y^{2} = I, \; Z^{2} = I
\end{align}
```

パウリ行列はエルミート、かつ、ユニタリです。


### X,Y,Z軸周りの回転行列

```math
\begin{align}
Rx(\theta) &= e^{i \theta X/2} =
	   \begin{pmatrix}
	   cos \frac{\theta}{2} & -i sin \frac{\theta}{2} \\
	   -i sin \frac{\theta}{2} & cos \frac{\theta}{2}
	   \end{pmatrix} \\
Ry(\theta) &= e^{i \theta Y/2} =
	   \begin{pmatrix}
	   cos \frac{\theta}{2} & -sin \frac{\theta}{2} \\
	   sin \frac{\theta}{2} & cos \frac{\theta}{2}
	   \end{pmatrix} \\
Rz(\theta) &= e^{i \theta Z/2} =
	   \begin{pmatrix}
	   e^{-i\theta/2} & 0 \\
	   0 & e^{i\theta/2}
	   \end{pmatrix} \\
\end{align}
```

行列表現がなぜこうなるかは、ユニタリ行列に対するオイラーの公式から、すぐにわかると思います。

### アダマール行列

定義は、以下です。

```math
H = \frac{1}{\sqrt{2}}
  \begin{pmatrix}
  1 & 1 \\
  1 & -1
  \end{pmatrix}
```

```math
H^{\dagger} = H, \; H^{2} = I
```

パウリ行列と同様、エルミート、かつ、ユニタリです。


### パウリ行列とアダマール行列に関する公式

```math
HZH = X, \; HXH = Z, \; HYH = -Y
```

### 回転行列とアダマール行列に関する公式

```math
HRz(\theta)H = Rx(\theta), \; HRx(\theta)H = Rz(\theta), \; HRy(\theta)H = Ry(-\theta)
```

### パウリ行列と回転行列に関する公式

```math
Rx(\pi/2)ZRx(-\pi/2) = Y, \; Z = Rx(-\pi/1)YRx(\pi/2)
```

また、ユニタリ行列に対するオイラーの公式を使うと、

```math
Rx(\pi/2)Rz(\theta)Rx(-\pi/2) = Ry(\theta)
```

を導くことができます。


### リー・トロッターの積公式

高校数学で習った、

```math
e^{A+B} = e^{A} \cdot e^{B}
```

は、$A,B$が行列で非可換($AB \neq BA$)の場合、成り立ちません。代わりに、十分大きい$N$を使って、以下の「リー・トロッターの積公式」で、近似的に表現できます。

```math
e^{A+B} \simeq (e^{\frac{A}{N}} \cdot e^{\frac{B}{N}})^{N}
```

だいたい、こんな前準備をしたところで、本題にうつります。


## 簡単な例

いきなり、一般的な数式を相手にするのはしんどいので、まず簡単なハミルトニアンの例で、量子回路を表現しながら、イメージをつかみたいと思います。

### 1粒子系

#### H = X

時間発展演算子は、

```math
U(t) = exp(-i X t) = Rx(-2t)
```

ですので、X軸周りの回転ゲート１個で表現できます。後の議論をすっきりさせるため、公式$Rx(\theta) = HRz(\theta)H$を使って、以下のようにゲート表現しておきます。

```
0 --H--Rz(-2t)--H--
```

#### H = Y

時間発展演算子は、

```math
U(t) = exp(-i Y t) = Ry(-2t)
```

ですので、Y軸周りの回転ゲート１個で表現できます。こちらも後の議論をすっきりさせるため、公式$Rx(\pi/2)Rz(\theta)Rx(-\pi/2) = Ry(\theta)$を使って、以下のようにゲート表現しておきます。

```
0 --Rx(-pi/2)--Rz(-2t)--Rx(pi/2)--
```

#### H = Z

時間発展演算子は、

```math
U(t) = exp(-i Z t) = Rz(-2t)
```

ですので、Z軸周りの回転ゲート１個で表現できます。

```
0 --Rz(-2t)--
```

### ２粒子系

簡単そうなパウリZからやってみます。

#### H = ZとZの直積

時間発展演算子は、

```math
U(t) = exp(-it Z_{0} \otimes Z_{1})
```

です。２粒子系(=２量子ビット)の状態は、$\ket{00},\ket{01},\ket{10},\ket{11}$なので、各々に$U(t)$を演算させてみます。

```math
\begin{align}
U(t) &= exp(-it Z_{0} \otimes Z_{1}) \ket{00} = exp(-it (1 \times 1)) \ket{00} = e^{-it} \ket{00} \\
U(t) &= exp(-it Z_{0} \otimes Z_{1}) \ket{01} = exp(-it (1 \times (-1))) \ket{01} = e^{it} \ket{01} \\
U(t) &= exp(-it Z_{0} \otimes Z_{1}) \ket{10} = exp(-it ((-1) \times 1)) \ket{01} = e^{-it} \ket{01} \\
U(t) &= exp(-it Z_{0} \otimes Z_{1}) \ket{11} = exp(-it ((-1) \times (-1))) \ket{11} = e^{it} \ket{11}
\end{align}
```

ちょっと考えれば、わかるのですが、以下のようなゲートを用意すれば、このような演算結果を再現することができます。

```
0 --*------------*---
1 --CX--Rz(-2t)--CX--
```

#### H = XとXの直積

時間発展演算子は、

```math
U(t) = exp(-it X_{0} \otimes X_{1})
```

です。$HZH = X$なので、

```math
\begin{align}
U(t) &= exp(-it X_{0} \otimes X_{1}) \\
&= exp(-it H_{0} Z_{0} H_{0} \otimes H_{1} Z_{1} H_{1}) \\
&= \cos(t)I - i \sin(t) (H_{0} Z_{0} H_{0} \otimes H_{1} Z_{1} H_{1}) \\
&= (H_{0} \otimes H_{1})(\cos(t)I - i \sin(t) Z_{0} \otimes Z_{1})(H_{0} \otimes H_{1}) \\
&= (H_{0} \otimes H_{1}) exp(-it Z_{0} \otimes Z_{1}) (H_{0} \otimes H_{1})
\end{align}
```

となります。最後の行の$exp$の部分は$H=Z_{0} \otimes Z_{1}$に対する時間発展なので、全体をゲートで表現すると、

```
0 --H--*------------*---H--
1 --H--CX--Rz(-2t)--CX--H--
```
となります。

#### H = YとYの直積

時間発展演算子は、

```math
U(t) = exp(-it Y_{0} \otimes Y_{1})
```

です。公式$Rx(\pi/2)ZRx(-\pi/2) = Y$を使うと、

```math
\begin{align}
U(t) &= exp(-it Y_{0} \otimes Y_{1}) \\
&= exp(-it Rx_{0}(\pi/2) Z_{0} Rx_{0}(-\pi/2) \otimes Rx_{1}(\pi/2) Z_{1} Rx_{1}(-\pi/2)) \\
&= Rx_{0}(\pi/2) \otimes Rx_{1}(\pi/2) exp(-it Z_{0} \otimes Z_{1}) Rx_{0}(-\pi/2) \otimes Rx_{1}(-\pi/2)
\end{align}
```

となります(省略しましたが途中でオイラーの公式を使っています)。最後の行の$exp$の部分は$H=Z_{0} \otimes Z_{1}$に対する時間発展なので、全体をゲートで表現すると、

```
0 --Rx(-pi/2)--*------------*---Rx(pi/2)--
1 --Rx(-pi/2)--CX--Rz(-2t)--CX--Rx(pi/2)--
```

となります。

#### H = XとYの直積

時間発展演算子は、

```math
U(t) = exp(-it X_{0} \otimes Y_{1})
```

です。これまでの議論と同様に、パウリ行列に関する公式とオイラーの公式を使って、

```math
\begin{align}
U(t) &= exp(-it X_{0} \otimes Y_{1}) \\
&= exp(-it H_{0} Z_{0} H_{0} \otimes Rx_{1}(\pi/2) Z_{1} Rx_{1}(-\pi/2)) \\
&= H_{0} \otimes Rx_{1}(\pi/2) exp(-it Z_{0} \otimes Z_{1}) H_{0} \otimes Rx_{1}(-\pi/2)
\end{align}
```

となります。最後の$exp$の部分は$H=Z_{0} \otimes Z_{1}$に対する時間発展なので、全体をゲートで表現すると、

```
0 --H----------*------------*---H---------
1 --Rx(-pi/2)--CX--Rz(-2t)--CX--Rx(pi/2)--
```

となります。

だいたいの傾向がわかってきたので、YとZの直積、ZとXの直積については、いきなりゲートで表現します。

#### H = YとZの直積

```
0 --Rx(-pi/2)--*------------*---Rx(pi/2)--
1 --I----------CX--Rz(-2t)--CX--I---------
```

#### H = ZとXの直積

```
0 --I--*------------*---I--
1 --H--CX--Rz(-2t)--CX--H--
```

となります。

### 多粒子系

まず、ハミルトニアンが、以下のようにパウリ演算子の直積のみで表わせる場合について考えます。

```math
H = P_{m} = \sigma^{p_{m}^{0}}_{0} \otimes \sigma^{p_{m}^{1}}_{1} \cdots \otimes \sigma^{p_{m}^{N-1}}_{N-1}

```

これまでの議論を地道にきちんとやっていくと、この時間発展は、以下のような量子ゲートとして表現できます。

```
0   --G--*-------------------------------*---G+--
1   --G--CX--*-----------------------*---CX--G+--
2   --G------CX--*----------------*--CX------G+--
     ...
N-1 --G------------CX--Rz(-2t)--CX-----------G+--
```

ここで、$G$は、対象となっている粒子に対するゲートが何かによって、切り替わります。

* $X$の場合、$G = H, \; G^{\dagger} = H$
* $Y$の場合、$G = Ry(-\pi/2), \; G^{\dagger} = Ry(\pi/2)$
* $Z$の場合、$G = I, \; G^{\dagger} = I$

です。

最後に、ハミルトニアンが、パウリ行列の直積を複数個足し合わせた形の場合です。これが最終的に求めたかったものになります。もう一度、ハミルとニアンの式を掲載します。

```math
\begin{align}
&H = \sum_{m=0}^{M-1} P_{m} \\
&P_{m} = \sigma^{p_{m}^{0}}_{0} \otimes \sigma^{p_{m}^{1}}_{1} \cdots \otimes \sigma^{p_{m}^{N-1}}_{N-1}
\end{align}
```

これに対する時間発展演算子は、

```math
U(t) = exp(-it \sum_{i=0}^{M-1} P_{m})
```

です。リー・トロッターの積公式を使うと、

```math
U(t) \simeq \bigl[\prod_{m=0}^{M-1} exp(-it P_{m}/D)\bigr]^{D}
```

と表せます(ここで、$D$は$t$に比べ十分大きな数です)。結局、どのような量子回路になるかというと、先程の$P_{m}$に対する量子回路を($t$を$t/D$に置き換えて)$P_{m}$の分だけ実行し、その全体を$D$回繰り返せば良いということになります。

## おわりに

定義や公式に基づき、なるべく自分で手を動かしながら計算してみました。間違っているところがあれば、ご指摘いただけるとうれしいです。





以上