[QuTiP]演算子の定義方法

# QuTiP?

pythonの計算ライブラリで，第二量子化したハミルトニアンの計算がしやすいものです。
[公式ページ](http://qutip.org/index.html)を見てみると，RIKENや学術振興会からサポートを受けているます。権威の傘を着るのはあまり好きではないですが，信頼性がありそうですね。
そもそも，このライブラリを使おうという研究者の方は英語が読めそうなので，ドキュメントの一部を日本語化した記事に需要があるかどうかは謎ですが，自分の備忘録として記載しておきます。

# インストール

[公式ドキュメント](http://qutip.org/docs/latest/installation.html)を読んでいると，`pip`と`conda`両方でのインストールは記述されていました。

`pip`であれば，

```
pip install qutip
```

`conda`であれば，

```
conda install -c conda-forge qutip
```
でインストールできます。

#　使い方

ここでの使い方は[公式ドキュメントの資料](https://nbviewer.jupyter.org/github/jrjohansson/qutip-lectures/blob/master/Lecture-0-Introduction-to-QuTiP.ipynb)のほとんどコピペです。

## 演算子の定義方法
QuTiPに計算させるObjectの作り方をまず説明します。

``` python3
from qutip import Qobj
```
Qobjというクラスを使って計算対象のオブジェクトを作ります。
これを使うことで適当なケット(列)ベクトルを生成できます。

``` python3
q = Qobj([[1], [0]])
```
と入力すると，jupyter notebook上では以下のような出力が得られます。

![スクリーンショット 2019-04-07 15.36.59.png](https://qiita-image-store.s3.amazonaws.com/0/195864/ad3c4510-2ceb-36f9-b9a9-66fecc73a9ad.png)

もちろんブラ（行）ベクトルに変更するのも簡単で，

```python3
q.dag()
```
と入力すると，出力が
![スクリーンショット 2019-04-07 15.42.16.png](https://qiita-image-store.s3.amazonaws.com/0/195864/21075fe4-f4e0-6299-5b16-7f9a03180c74.png)
と得られます。つまり，ケット→ブラベクトルへの変換を数式から直感的に操作可能です。

試しに，

```python3
q.dag()*q
q*q.dag()
```
の計算を行ってみると，それぞれ

![スクリーンショット 2019-04-07 15.44.18.png](https://qiita-image-store.s3.amazonaws.com/0/195864/bacedca2-05ab-b1c6-cd82-2ada7c1794b7.png)
![スクリーンショット 2019-04-07 15.44.23.png](https://qiita-image-store.s3.amazonaws.com/0/195864/6d3cb41a-7f40-46d4-6376-25310268b89f.png)
と計算できるので，簡単な記述でこれらの計算を行ってくれることがわかります。

さて，スピン2状態のみを考えたときにそのエネルギーを表すハミルトニアンは，

```math
\mathcal{H}=\frac{\hbar\omega}{2}\sigma_z=\frac{\hbar\omega}{2}\{|+\rangle\langle +| - |-\rangle\langle -|\}
```
と記載されます。ここで，$|+\rangle$はエネルギーが高い方の状態，$|-\rangle$はエネルギーが低い方の状態，と定義します。また，先程定義した`q`は$|+\rangle$と一致します。
では，このハミルトニアンの演算子の部分($\sigma_z$)をQuTiPを使って定義してみましょう。
Qobjクラスは行列にも対応しており，以下のように定義することで，この演算子を定義することができます。
```python3
sigma_z = Qobj([[1,0],[0,-1]])
```
![スクリーンショット 2019-04-07 16.00.26.png](https://qiita-image-store.s3.amazonaws.com/0/195864/3adb1e34-7259-4a2d-eed4-0d3914c80668.png)

なので，あとはハミルトニアンを作るのであれば適当な係数をつければ良いだけなので簡単ですね。

他の系との相互作用を考えるのであれば，$\sigma_y$や$\sigma_x$が必要です。
複素数も簡単に定義できるので，$\sigma_y$を定義してみます。

```python3
sigma_y = Qobj([[0,-1j],[1j,0]])
```
![スクリーンショット 2019-04-07 16.10.27.png](https://qiita-image-store.s3.amazonaws.com/0/195864/42ca2739-85dd-f119-c2db-9cc320976b3f.png)
このように，`-1j`と書くことで複素数の定義も簡単に記述できます。

これら定義した演算子の次元は

```python3
hoge.dims
```
で確認できます。例えば，上で定義した$\sigma_y$に対してやってみましょう。

```python3
sigma_y.dims
```

```output
[[2], [2]]
```

## 内部に定義されているクラスを使う
今までほQobjクラスをそのまま使っていましたが，よく利用される演算子についてはQuTiP内で定義されいることがあります。これらを使うとより簡単な記述で定義ができます。

### 二準位系の演算子

二準位系の状態ベクトルを記載するだけであれば，専用のクラスがあります。
例えば，$\sigma_z$と$\sigma_y$であれば，

```python3
from qutip import sigmay, sigmaz
sigma_y = sigmay()
sigma_z = sigmaz()
```
で上記と同じものが定義できます。

### Fock stateの状態ベクトル
Fock stateという数状態の状態ベクトルを作ります。
Fock stateは例えばフォトンの数で基底を張る状態です。フォトンが1個であれば$|1\rangle$，10個であれば$|10\rangle$です。

```python3
from qutip import basis, fock
basis(4,2)
fock(4,2)
```
この２つであれば今はあまり差はありません。
両方同じ出力が得られます。

![スクリーンショット 2019-04-07 16.24.05.png](https://qiita-image-store.s3.amazonaws.com/0/195864/eaf6b3f3-1231-4c4a-aea3-83da242e8dd7.png)

`qobj`クラスを使って頑張って定義してもよいですが，こちらのほうが簡単に定義できますね。

Fock状態の生成消滅演算子の行列もサポートされています。

```python3
from qutip import destroy,create
N = 4
a = destroy(N)
a_dag = create(N)
```
ここで`N`はFock状態の次元です。ヒルベルト空間上にどこまで基底を張るのかをここで定義しています。
さて，定義した`a`を出力してみると，
![スクリーンショット 2019-04-07 16.33.08.png](https://qiita-image-store.s3.amazonaws.com/0/195864/81dabe59-17ce-670f-dc4d-4e79b476d17a.png)
であるので，よくある文献と一致しています。便利ですね。

せっかくなので先程記述したFock状態ベクトルにここで定義した生成，消滅演算子を作用させてみましょう。

#### 生成演算子
```python3
adag*fock(4,2)
```
![スクリーンショット 2019-04-07 16.35.06.png](https://qiita-image-store.s3.amazonaws.com/0/195864/5b731451-e3ef-f222-7c0f-db0059a2a9e5.png)
生成演算子の定義から

```math
a^\dagger|n\rangle = \sqrt{n+1}|n+1\rangle
```

であるので正しく計算ができていることがわかります。

#### 消滅演算子

```python3
a*fock(4,2)
```
![スクリーンショット 2019-04-07 16.38.14.png](https://qiita-image-store.s3.amazonaws.com/0/195864/94bf3fb5-fcec-c566-f6ee-ee596b9430a5.png)

これも消滅演算子の定義から

```math
a|n\rangle = \sqrt{n}|n-1\rangle
```

であるため，正しく計算できていることがわかります。

# まとめ

今回はQuTiPでの計算を行うために，演算子の定義をいくつか行いました。
[公式ドキュメント](https://nbviewer.jupyter.org/github/jrjohansson/qutip-lectures/blob/master/Lecture-0-Introduction-to-QuTiP.ipynb)により詳細は書かれているので，そちらを参照してください。
次は相互作用がある場合を考えて，もうちょっと面白い計算を行います。