QuTIPとipywidgetsを組み合わせてみるとモデルの理解に役立つ、という話

# 緒言

以前、QuTIPという量子ダイナミクスの計算のためのフレームワークについて、簡単な[解説記事](https://qiita.com/swakamoto/items/1226b0d53af8b287b833)を書いた。
QuTIPの[チュートリアル](http://qutip.org/tutorials.html)には、これ以外にも多くの例が挙げられており、ライブラリの使い方の理解などに大変役立つ。

これらは単純にcloneして動かして見ても良いのだが、ipywidgetsというライブラリを使うことでjupyter notebook上でパラメータなどをインタラクティブに操作しながら結果のグラフを表示させたりすることができ、モデルの理解に有用だと思う。
本記事では、これについて紹介したい。

# ipywidgetsのインストールなど
ipywidgetsのインストールは、pipで簡単に行うことができる。
解説などは、すでに有用なページ[^1]などがあるので、それらを参考にしてほしい。

# Jammis-Cummingsモデルをinteractiveにしてみる
チュートリアルに紹介されている[Jaynes-Cummingsモデル](https://nbviewer.jupyter.org/github/qutip/qutip-notebooks/blob/master/examples/rabi-oscillations.ipynb)を動かしてみる。
これをipywidgetsを使って書き直して見たコードは、[ここ](https://nbviewer.jupyter.org/github/YukiSakamoto/QuTIP_models/blob/master/Interactive_Rabi_Oscillation.ipynb)にある。

やることは簡単で、チュートリアルではトップレベルに書かれているシミュレーションのプログラムを関数化して、それを`interact()`関数で呼んでやるだけだ。

```python
def calculate(wc_per_2pi, wa_per_2pi, g_per_2pi, kappa, gamma, n_th_a, use_rwa = True):
    wc = wc_per_2pi  * 2 * np.pi  # cavity frequency
    wa = wa_per_2pi  * 2 * np.pi  # atom frequency
    g  = g_per_2pi * 2 * np.pi  # coupling strength
    #kappa = 0.005          # cavity dissipation rate
    #gamma = 0.05           # atom dissipation rate
    N = 15                 # number of cavity fock states
    #n_th_a = 0.0           # temperature in frequency units
    #use_rwa = True

    tlist = np.linspace(0,25,100)
    
    # intial state
    psi0 = tensor(basis(N,0), basis(2,1))    # start with an excited atom

    # operators
    a  = tensor(destroy(N), qeye(2))
    sm = tensor(qeye(N), destroy(2))

    # Hamiltonian
    if use_rwa:
        H = wc * a.dag() * a + wa * sm.dag() * sm + g * (a.dag() * sm + a * sm.dag())
    else:
        H = wc * a.dag() * a + wa * sm.dag() * sm + g * (a.dag() + a) * (sm + sm.dag())
    
    c_op_list = []

    rate = kappa * (1 + n_th_a)
    if rate > 0.0:
        c_op_list.append(np.sqrt(rate) * a)

    rate = kappa * n_th_a
    if rate > 0.0:
        c_op_list.append(np.sqrt(rate) * a.dag())

    rate = gamma
    if rate > 0.0:
        c_op_list.append(np.sqrt(rate) * sm)
    
    output = mesolve(H, psi0, tlist, c_op_list, [a.dag() * a, sm.dag() * sm])
    
    fig, ax = plt.subplots(figsize=(8,5))
    ax.plot(tlist, output.expect[0], label="Cavity")
    ax.plot(tlist, output.expect[1], label="Atom excited state")
    ax.legend()
    ax.set_xlabel('Time')
    ax.set_ylabel('Occupation probability')
    ax.set_title('Vacuum Rabi oscillations');

from ipywidgets import interact
interact(calculate, wc_per_2pi = (0, 2, 0.1), wa_per_2pi = (0, 2, 0.1), g_per_2pi = (0.00, 0.1, 0.01), 
         kappa = (0.0, 0.01, 0.001), gamma = (0.0, 0.10, 0.01), n_th_a = 0 )

```

# 動かしてみよう
jupyter notebookでこれを実行すると、以下のような画面が現れる。
パラメータを変化させるためのバーが出ている。なお、この初期設定パラメータは全て公式のチュートリアルにあるパラメータと同じ値である。

![スクリーンショット 2019-06-10 22.00.48.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/24675/44193b0f-a840-1035-93d2-76932d52f6c4.png)

これをグリグリやってみる。
<blockquote class="twitter-tweet" data-lang="ja"><p lang="ja" dir="ltr">QuTIPとJupyter-notebookで遊んで見た <a href="https://t.co/bUgG6cnKQx">pic.twitter.com/bUgG6cnKQx</a></p>&mdash; swakamoto (@swakamoto) <a href="https://twitter.com/swakamoto/status/1138072284348813313?ref_src=twsrc%5Etfw">2019年6月10日</a></blockquote>
<script async src="https://platform.twitter.com/widgets.js" charset="utf-8"></script>


こうやってみると、「共振器と原子の振動数(wcとwa)が近接していないと、振動しないのだな。。。」とか、「gの値（共振器と原子の結合定数）が大きくなると、振動数が大きくなる」などといったことがわかりやすく、自分で引数をぽちぽち変えるよりもモデルの性格などが見えて来やすいと思う。

この手の数値シミュレーションにおいては、QuTIPに限らなくても有用な場面は多いと思うので、今後も「まずは遊んでみるつもりで」このツールを積極的に使っていきたいと思った次第である。

[^1]: たとえば、https://kivantium.net/python-ipywidgets など