🪐

CUDA-Q で遊んでみる (3) — GPU で 期待値計算をしてみる

に公開

目的

CUDA-Q で遊んでみる (2) — GPU で Grover のアルゴリズムを回してみる に続いて別の用途で GPU を使ってみようというもの。

面白いアプリケーションを思いつかなかったので、非常にくだらない計算をする。今度は、cuQuantumテンソルネットワークシミュレータの cuTensorNet を間接的に使用する。

以前にも 行列積状態について考える (10) — 50 量子ビットの期待値計算 で「行列積状態」と呼ばれるテンソルネットワークの一種を用いて大量の量子ビットでの期待値計算を行った。cuTensorNet はこの種の計算に向いているので今回確認してみたい。最終的に 200 量子ビットの期待値計算を行う。

量子回路の実装をする

大して面白い例が思いつかなかったことから、さくっと実装して済ませる。

まずは必要なパッケージをインストールする

今回は CUDA-Q だけを使う。

!pip install -qU cudaq"

アスキーアートの回路図がちゃんと表示されるようにする

ユーティリティを定義する。

from IPython.display import HTML, display

def show_fixed(text, font="Consolas, Roboto Mono, monospace", size=13):
    text = text.expandtabs(4)   # タブをスペースに変換
    esc = (text.replace("&", "&")
               .replace("<", "&lt;")
               .replace(">", "&gt;"))
    html = f'<pre style="font-family:{font}; font-size:{size}px; white-space:pre; font-variant-ligatures:none;">{esc}</pre>'
    display(HTML(html))

必要なパッケージをインポートする

今回は期待値を計算するハミルトニアンを作るので spin もインポートする。

import numpy as np
import cudaq
from cudaq import spin

CUDA-Q カーネルを定義する

@cudaq.kernel
def kernel(num_qubits: int):
    qubits = cudaq.qvector(num_qubits)
    x(qubits[0])
    for i in range(1, num_qubits):
        x.ctrl(qubits[0], qubits[i])

try:
    kernel.compile()
except Exception as e:
    print(e)

可視化する

何かをしたかったのだけど、常に期待値が 0 になる回路になってしまって、無理やり調整した結果、本当に馬鹿げた回路になった。

try:
    show_fixed(str(cudaq.draw(kernel, 5)))
except Exception as e:
    print(e)

期待値を計算する

以下のようなハミルトニアンの係数を作成する。

np.random.seed(42)

# 20: 15.7 s
# 21: 34.8 s
# 22: 1min 17s
num_qubits = 21

coefficients = np.random.random(num_qubits).tolist()

np.sum(coefficients).item()

9.765898662678088

上記の係数を使ってハミルトニアンを以下のように定める。前述の量子回路のもとでのこのハミルトニアンの期待値は、係数の和「9.765898662678088」にマイナスをつけたものになる。

operator = sum([coefficients[i] * spin.z(i) for i in range(num_qubits)])

print(operator)

(0+0i) + (0.37454+0i) * Z0 + (0.950714+0i) * Z1 + (0.731994+0i) * Z2 + (0.598658+0i) * Z3 + (0.156019+0i) * Z4 + (0.155995+0i) * Z5 + (0.0580836+0i) * Z6 + (0.866176+0i) * Z7 + (0.601115+0i) * Z8 + (0.708073+0i) * Z9 + (0.0205845+0i) * Z10 + (0.96991+0i) * Z11 + (0.832443+0i) * Z12 + (0.212339+0i) * Z13 + (0.181825+0i) * Z14 + (0.183405+0i) * Z15 + (0.304242+0i) * Z16 + (0.524756+0i) * Z17 + (0.431945+0i) * Z18 + (0.291229+0i) * Z19 + (0.611853+0i) * Z20

実際に計算を実行すると、丸め誤差を除いて宣言通りに「9.765898662678088」にマイナスをつけたものが得られる。

%%time

try:
    result = cudaq.observe(kernel, operator, num_qubits)
    expval: float = result.expectation()
    print(f"{expval=}")
except Exception as e:
    print(e)

expval=-9.76589866267809
CPU times: user 33.8 s, sys: 501 ms, total: 34.3 s
Wall time: 34.8 s

ネタバレとしては、以下を計算したのである。

\begin{align*} \braket{1^{\otimes n}|H|1^{\otimes n}} &= \braket{1^{\otimes n}|\sum_{i=1}^n \alpha_i Z_i|1^{\otimes n}} \\ &= \sum_{i=1}^n \alpha_i \braket{1^{\otimes n}|Z_i|1^{\otimes n}} = \sum_{i=1}^n \alpha_i \braket{1|Z|1} = -\sum_{i=1}^n \alpha_i \end{align*}

ところでこの計算は遅い。コメントとしてコードに埋め込んでいるが、20 量子ビットで 15.7 s、21 量子ビットで 34.8 s、22 量子ビットで 1min 17s くらいの遅さである。

そこで、

GPU シミュレーションを行う

ここからは GPU (cuTensorNet) の出番である。NVIDIA T4 以上が使用可能な Colab で起動している ものとし、どこまで高速化できるかを確認する。

ターゲットデバイスを GPU(テンソルネットワーク)にセットする

追加の設定は以下を実験の前に実行するだけである。

if cudaq.num_available_gpus() > 0:
    cudaq.set_target('tensornet')

他の実装は上で書いた通りをそのまま流用できるので、早速実験してみる。

np.random.seed(42)

# 20: 432 ms
# 21: 460 ms
# 22: 510 ms
# 30: 943 ms
# 40: 1.92 s
# 50: 2.42 s
# 60: 3 s
# 70: 4.11 s
# 80: 5.25 s
# 90: 9.46 s
# 100: 11 s
# 125: 15.9 s
# 150: 21.3 s
# 175: 32.4 s
# 200: 45.8 s
num_qubits = 21

coefficients = np.random.random(num_qubits).tolist()

np.sum(coefficients).item()

9.765898662678088

operator = sum([coefficients[i] * spin.z(i) for i in range(num_qubits)])

print(operator)

(0+0i) + (0.37454+0i) * Z0 + (0.950714+0i) * Z1 + (0.731994+0i) * Z2 + (0.598658+0i) * Z3 + (0.156019+0i) * Z4 + (0.155995+0i) * Z5 + (0.0580836+0i) * Z6 + (0.866176+0i) * Z7 + (0.601115+0i) * Z8 + (0.708073+0i) * Z9 + (0.0205845+0i) * Z10 + (0.96991+0i) * Z11 + (0.832443+0i) * Z12 + (0.212339+0i) * Z13 + (0.181825+0i) * Z14 + (0.183405+0i) * Z15 + (0.304242+0i) * Z16 + (0.524756+0i) * Z17 + (0.431945+0i) * Z18 + (0.291229+0i) * Z19 + (0.611853+0i) * Z20

準備が出来たので期待値を計算する

期待値を計算する

%%time

try:
    result = cudaq.observe(kernel, operator, num_qubits)
    expval: float = result.expectation()
    print(f"{expval=}")
except Exception as e:
    print(e)

expval=-9.76589866267809
CPU times: user 458 ms, sys: 2.99 ms, total: 461 ms
Wall time: 460 ms

驚異的に速くなってしまった。CPU 計算で 34.8 s だったものが、GPU を活用したテンソルネットワーク計算で 460 ms になってしまった。

量子ビット数と経過時間の関係を見る

折角なので量子ビット数を色々と変更して経過時間を観察してみた。200 量子ビットでも僅か 46 秒程度 である。勿論量子回路の内容やハミルトニアンにもよると思うが、1 つの目安にはなるであろう。

こういう感じのグラフはとりあえず対数グラフにして見たい。何とも言えないが、量子ビット数が十分に大きいところ (\geq 90) では直線状になっていそうな気もするので、経過時間は指数的に増加しそうな気がする。

まとめ

CUDA-Q を用いた期待値計算をやってみた。内容は大変つまらないものになったが、爆速で計算できそうな感触は得た。

但し、量子回路が表現しているシステムのエンタングルメントエントロピー (Entropy of entanglement) が十分に小さい場合でないと爆速さは満喫できない。何故ならエンタングルメントが増大するとノード間のエッジの結合次元が増大するので、いわゆる縮約計算のコストが増大するからである。このため、実験をしたいシステムのエントロピーや結合次元が、量子ビット数に対してどのようにスケールするかを事前に確認しておくことは重要であろう。場合によっては CUDA-Q で遊んでみる (2) — GPU で Grover のアルゴリズムを回してみる で触れた cuStateVec のほうが計算に適している可能性もある。

次は QAOA にでもチャレンジしてみたい。

参考文献

GitHubで編集を提案

Discussion