量子計算のQmod言語の仕様を眺めて少し触ってみた感想

こんにちは! Insight Edge分析チームの梶原(悠)です。最近ひょんな経緯で量子計算用のQmodという言語のフィジビリ兼ゆる勉強会に顔を出しています。 Qmod言語 1 はclassiq社という量子ベンチャーが提供している無償有償ツール※で、簡便に量子アルゴリズムを実装できる高水準言語をうたっています。私は量子計算について何も知らない素人ですが、基本的なpythonと線形代数の知識があれば使えるとのことで、量子畑の人たちにあれこれ教えていただきながら、すこし触ってみました。量子計算に興味や前提知識はないが、技術動向はある程度把握しておきたいと考える技術者の読み手を想定して、言語仕様の一部やツールに触れてみた感想などを書きます。
※ 2025.3.7 classiq社様より無償ツールの記述は誤りであるとのご連絡を頂いた為訂正しました。公式サイト中の"free for non-comercial purpose"の記載に基づき無償ツールと解釈していたのですが、最近有償化されたそうです。この誤りを含む複数のご指摘に感謝申し上げます。

はじめに

量子計算とは

古典的なコンピュータにおけるCPUのレジスタは一刻にひとつの状態しか取れません。これは重ね合わせ状態の崩れた物理系を利用して実装されているためです。一方、重ね合わせ状態にある物理系でレジスタを実装すれば、指数的に多くの状態を並列して時間発展させることができます。量子計算では、指数的に大きな並列処理を利用して、複雑な計算を高速に実行します。ただし結果の観測などに強い制約があり、どんな問題でも指数的に速く解けるわけではありません。

実用化に向けた現状

ibmroadmap

IBM社のロードマップ 2 ではエラーの影響なく量子計算を実行できるデバイスの提供は2029年以降とされています。 classiq社のセミナー講師の方の話では、量子計算には実ビジネスで意味のある活用例はまだ存在しないようです。何かのイベントの質疑応答で理研の研究者の方も「量子乱数は実用されているが、最適化などの文脈で意味のある活用例はまだない」と回答されてました。これらの話を聞くと、量子技術を付加価値にして実社会に役立つものづくりをすることは現時点で難しいように思えます。

しかし、興味深いベンチャーも出てきています。量子計算が実用化すると楕円曲線暗号を現実的な時間内に解けるようになるため、既存の暗号通貨が攻撃される懸念があります。 BlocQ 3 というベンチャー企業では、暗号生成部分のみ量子計算を利用したブロックチェーン技術を採用し、量子技術の実用化が進んでも進まなくてもワークする量子readyな暗号通貨を開発しているようです。量子に限らず、まだ早すぎる技術をどう売るか考える際に、〇〇readyな〜という切り口は学ぶ処がありそうに感じました。

量子回路モデル

quantumcomputing

Qmodは量子回路という計算モデルの記述言語です。量子回路を提案したD.Deutschは計算器を次のように説明しています。 4

直感的には, 計算器とは, 系の動的な発展で入力状態の集合から出力状態の集合へと遷移する任意の物理系です. 各状態には, なんらかの自然なラベルがつけられていて, マシンは指定された入力ラベルの状態で準備されます. そして, いくらかの動作を経たのち, 出力状態のラベルが観測されます.

古典計算では、入出力のラベルやレジスタの状態は01列で表されていました。量子計算でもラベルは01列ですが、レジスタは複素ベクトルで表されます。古典計算の時間発展はレジスタの状態をコピーや論理演算子などでくりかえし変換して記述されました。量子計算はレジスタにゲートと呼ばれるユニタリ行列をくりかえし作用させて時間発展を記述します。時間発展は可逆でなければならず、状態のコピーはできません。レジスタに作用させるゲートを並べた系列を量子回路 5 といいます。下記の表に量子回路による計算の構成要素をまとめます。

用語	説明
キュビット (qubit)	大きさが1の2次元複素空間のベクトル $\langle \psi \| \in \mathbb{C}^2$ をキュビットという.
量子レジスタ	大きさが1のn次テンソル空間のベクトル $\langle \psi \| \in (\mathbb{C}^2)^{\otimes n}$ を長さnの量子レジスタという.
キュビットの状態 6	トレースが1の半正定値なエルミート行列を密度行列という. キュビット $\langle \psi \| \in \mathbb{C}^2$ で定まる密度行列 $\| \psi \rangle \langle \psi \|$ をキュビットの状態という.
量子レジスタの状態	量子レジスタ $\langle \psi \| \in (\mathbb{C}^2)^{\otimes n}$ で定まる密度行列 $\| \psi \rangle \langle \psi \|$ を量子レジスタの状態という.
量子ゲート	量子レジスタ上のユニタリ行列を量子ゲートという.
量子回路	長さ $𝑚$ のレジスタ上の量子ゲート $G$ は自然に長さ $n \geq 𝑚$ のレジスタ上のユニタリ行列 $U$ に拡張できる. 量子ゲート $G _ 0,...,G _ T$ を長さ $𝑛$ のレジスタ上に拡張したユニタリ行列の列 $U _ 0,…,U _ 𝑇$ を𝑛ビット量子回路とする.
観測可能量	冪等 (すなわち $P ^ 2=P$ ) な行列 $P$ を、部分空間 $Codom(P)$ への射影作用素という. エルミート行列 $\mathfrak{o} _ j = \mathbf{I}^{\otimes j -1} \otimes \| 1 \rangle \langle 1 \| \mathbf{I}^{\otimes n-j}$ を $𝑗$ 番キュビットのオブザーバブルという. オブザーバブル $\mathfrak{o} _ 𝑗$ の固有値 $𝑦 _ 𝑗 \in \{ 0,1 \}$ に対応する固有空間へのエルミートな射影作用素を $P _ j (y _ j)$ とかく.
測定	量子レジスタ $\langle \psi \| \in (\mathbb{C}^2)^{\otimes n}$ の $𝑗 _ 1,…,𝑗 _ 𝑛$ キュビットの状態の測定とは以下の1. 及び2. を行うことである. 1. (統計公式) ビット列のサンプリング $bits _ {out} \sim Pr(bits _ {out} \| bits _ {in})=tr(P \| \psi \rangle \langle \psi \|)$ を行う. ただしビット列 $bits _ {out}=(y _ {j _ 1},…,y _ {j _ m })$ に対応する射影作用素を $P=P _ {j _ 1} (y _ {j _ 1}) ... P _ {j _ 1} (y _ {j _ m})$ とした. 2. (射影仮説) レジスタの状態を密度行列 $\| \psi \rangle \langle \psi \|$ から $P \| \psi \rangle \langle \psi \| P / tr(P \| \psi \rangle \langle \psi \| P)$ に変える.

量子回路は下図のような直感的なダイアグラムで記述することができます。横線はキュビットを表し、四角い箱はゲートを表します。横線の左端がキュビットの初期状態に対応し、左から順にゲートを作用させていくことで時間発展を記述します。

quantumcircuit

Qmod言語の仕様調査

コンパイル階層

layerclassiq

量子計算を実行する物理デバイスによってキュビットの結合状況や利用できる基本ゲートの種類は異なります。 qmodツールは抽象的な量子回路を物理デバイス上で実行可能な量子回路へとコンパイルします。 qmodのコンパイル階層 7 の詳細は不明ですが、概ね上図のような流れと想像されます。

QNum型

qmodにはいくつかのQuantum型が定義されています。 QBit型はキュビットに対応します。 QNum型は量子レジスタに実数の固定少数点表現として解釈するためのメタ情報を付与したものと思われます。いくつかのパタンを実験した限りでは基底ベクトルと実数は下記のように対応するようです。

fixedpointfloat

量子レジスタの解放

qmodは量子レジスタの領域の割当や解放を自動的に行います。 Quantum型の変数に割り当てられた領域を明示的に解放させる方法として少なくとも以下の3通りはあるようです.

bind(src,dst)関数のsrc側の引数にする.
qfuncデコレータのついた関数のinput修飾子のついた引数にする.
within_apply(within,apply)関数のwithin部で宣言、初期化する。

ライブラリ関数の例: 標準ゲート関数(X, Y, Z)

X,Y,Z関数はそれぞれ以下で定義されるパウリ $X,Y,Z$ 演算子を表します。

#	分類・名称	記号	定義	基底 $\{ \|0 \rangle, \|1 \rangle \}$ についての行列表示
1	パウリX演算子, NOT	$X: \mathbb{C}^2 \rightarrow \mathbb{C}^2$	$X=\| 0 \rangle \langle 1 \| + \| 1 \rangle \langle 0 \|$	$X = \begin{pmatrix} 0 & 1 \\ 1 & 0 \end{pmatrix}$
2	パウリY演算子	$Y: \mathbb{C}^2 \rightarrow \mathbb{C}^2$	$Y=-i \| 0 \rangle \langle 1 \| + i \| 1 \rangle \langle 0 \|$	$Y = \begin{pmatrix} 0 & -i \\ i & 0 \end{pmatrix}$
3	パウリZ演算子	$Z: \mathbb{C}^2 \rightarrow \mathbb{C}^2$	$Z=\| 0 \rangle \langle 0 \| - \| 1 \rangle \langle 1 \|$	$Z = \begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix}$

ライブラリ関数の例: linear_pauli_rotations関数

linear_pauli_rotations関数は以下の形で呼び出します。

$\displaystyle{ \begin{align*} linear \_ pauli \_ rotations( \\ &bases = [\mu _ 1,...,\mu _ m ], \\ &slopes = [a _ 1,...,a _ m ],\\ &offsets = [b _ 1,...,b _ m ],\\ &x = |x \rangle _ n ,\\ &q = | q _ 1 q _ 2... q _ m \rangle _ m\\ ) \end{align*} }$

(ただし、 $\mu _ 1,...,\mu _ m \in \{ Pauli.X.value, Pauli.Y.value, Pauli.Z.value \}$ はパウリ行列を指定するenum。)

この関数は以下のようにレジスタ $|q _ 1 q _ 2 ... q _ m \rangle _ m$ を変化させます.

Input	Output
$\|x\rangle _ n \otimes \left( \bigotimes _ {k=1,...,m} \|q _ k \rangle \right)$	$\|x\rangle _ n \otimes \left( \bigotimes _ {k=1,...,m} \exp \left( -i \left( \frac{a _ k}{2} x + \frac{b _ k}{2} \right) P _ k \right) \|q _ k \rangle \right)$

(ただし $P _ k \in \{X, Y, Z \}$ は引数basesで指定したパウリ行列.)

linear_pauli_rotations関数はキュビットを行列 $\exp(i \theta P)$ で変換する関数だとわかりました。そこで行列 $\exp(i \theta P)$ による変換が何を意味するかを考えてみます。

hophmap

写像 $f: \mathbb{C} ^ 2 \rightarrow \mathbb{R} ^ 3$ を $f( | \psi \rangle) = \left( \langle \psi | X | \psi \rangle, \langle \psi | Y | \psi \rangle, \langle \psi | Z | \psi \rangle \right)$ で定め、ホップ写像と呼びます。ホップ写像を $\mathbb{C} ^ 2$ の単位球に制限すると $\mathbb{R} ^ 3$ の単位球への全射を与えます。この全射の像をブロッホ球と呼びます。キュビット $| \psi \rangle$ を行列 $\exp(i \theta P)$ で変換すると、ブロッホ球上の点 $f(| \psi \rangle)$ は点 $f(\exp(i \theta P) | \psi \rangle )$ へ動きます。ホップ写像の定義式にパウリ行列 $P \in \{X, Y, Z \}$ の指数を三角関数で表した $\exp(i \theta P) = \cos(\theta) \mathbf{I} + i \cdot \sin(\theta) P$ を代入すると次の関係が導かれます。

$\displaystyle{ f(\exp(i \theta X) \ | \psi \rangle )^t = \begin{pmatrix} 1 & 0 & 0 \\ 0 & \cos2\theta & -\sin2\theta \\ 0 & \sin2\theta & \cos2\theta \\ \end{pmatrix} f(| \psi \rangle )^t }$

したがって、行列 $\exp(i \theta X)$ はブロッホ球を $x$ 軸周りに $2\theta$ 回転させます。また、 $X,Y,Z$ をそれぞれ $Y,Z,X$ にぐるりと入れ替えて同様の議論を繰り返せば、行列 $\exp(i \theta Y)$ が $y$ 軸周りの $2\theta$ 回転を引き起こすことがわかります。同様にして、行列 $\exp(i \theta Z)$ は $z$ 軸周りの $2\theta$ 回転を引き起こすことがわかります。

回路生成の試行

量子モンテカルロという手法で積分 $\int _ {0} ^ {\frac{\pi}{6}}\int _ {0} ^ {\frac{\pi}{6}} \sin ^ 2(x+y)dxdy$ の値を推定する回路を題材に、回路生成を試行してみます。量子モンテカルロでは、積分値の推定問題をユニタリ行列の固有値問題に帰着して解きます。大まかなアイディアとしては、まず推定したい積分値の情報を固有値として持つユニタリ行列を構成します。次に位相キックバックというトリックを用いて固有値をキュビットの位相として取り出します。最後に、位相を量子フーリエ変換の逆変換等で振幅に反映させて観測すれば、興味のある積分値を高速に推定できるという流れです。推定したい積分値の情報を固有値として持つユニタリ行列を構成するには2つの材料が必要です。ひとつめは、求めたい積分値に応じた角度 $\theta$ だけ初期ベクトルを回転させながらある実2次元部分空間 $V$ に送り込む作用素(state_loading)です。ふたつめは、部分空間 $V$ を第一軸について鏡映する作用素(good_state_oracle)です。これらの材料を掛け合わせることでGroverOperatorという作用素を構成できます。 GroverOperatorは空間 $V$ を $2\theta$ 回転させるため、固有値が積分値の情報を持つ作用素を作れたことになります。実装上は、材料となる2つの作用素(state_loading及びgood_state_oracle)を qmodの組み込み関数amplitude_estimationに与えるだけで、 qpeによる量子モンテカルロの回路を作成できます。以下は上記の積分値を推定する回路を生成してtranspile前後の回路をOpenQASMという形式で取得する実装例です。

    from classiq import * # バージョン 0.67.0 で試行.

    n_qpe = 5
    D = 2
    n_prn = 5
    a = 2.0 * (math.pi / 6.0) / (2 ** n_prn)
    b = 0.0
    
    def good_state_oracle(state: QNum[QBit]):
        x = QNum("x_", n_prn*D)
        q = QBit("q_")
        bind(state, [q, x])
        Z(q)
        bind([q, x], state)
        return state

    def state_loading(state: QNum[QBit]):
        x1 = QNum("x1", n_prn)
        x2 = QNum("x2", n_prn)
        q = QBit("q")
        bind(state, [q, x1, x2])
        
        hadamard_transform(x1)
        linear_pauli_rotations(
            bases=[Pauli.Y.value], slopes=[a], offsets=[b], x=x1, q=q # sin((ax+b)/2)
        )
        
        hadamard_transform(x2)
        linear_pauli_rotations(
            bases=[Pauli.Y.value], slopes=[a], offsets=[b], x=x2, q=q # sin((ax+b)/2)
        )
        bind([q, x1, x2], state)
        return state

    @qfunc
    def main(phase: Output[QNum]):
        state = QNum[QBit]("state")
        allocate(n_prn*D+1, state)
        allocate_num(n_qpe, False, n_qpe, phase)
        amplitude_estimation(
            oracle=good_state_oracle, # S_psi
            space_transform=state_loading, # A
            phase=phase,
            packed_vars=state
        )

    model = create_model(main)
    model = set_constraints(
        model,
        Constraints(
            max_width=25, 
            optimization_parameter=OptimizationParameter.DEPTH
    ))
    
    model = set_preferences(model, Preferences(
        optimization_timeout_seconds = 3600,
        timeout_seconds = 3600 * 2,
        qasm3=True,
        custom_hardware_settings = CustomHardwareSettings( 
            basis_gates=["rz", "cx", "sx", "x"],
        ),
        transpilation_option="intensive",
        optimization_level = 3,
        debug_mode=False
    ))
     
    qprog = synthesize(model)
    circuit = QuantumProgram.from_qprog(qprog)
    qasm_transpiled = circuit.transpiled_circuit.qasm
    qasm_raw = str(circuit.qasm)

上記の回路生成を回路幅の制約や最適化有無を切り替えながら何度か実行してみました。 synthesize関数が生成したtranspile前後の回路の幅(キュビットの個数)と深さ(時間発展のステップ数)の関係のグラフを示します。最適化なしの設定ではすべての回路が同じ幅になりました。最適化ありの設定で生成したtranspile後の回路には幅と深さのトレードオフが見られます。一方、transpile前の回路で幅の制約を緩めていくと、深さはそのままで幅だけが広がっていくことが観察されました。これは敢えて幅を広げた冗長な回路にすることで、transpile後の回路がより圧縮されるように準備しているものと想像します。 transpile後の回路の深さは幅23ほどで改善が頭打ちになり、その後はゲート数だけが悪化しています。仕様上、幅と深さは同時に最適化できないため、最適な回路を生成するためには、ユーザ側で回路幅の適切な制約を探索する必要がありそうです。

depthcomparison

感想

今回は量子計算用のqmod言語の仕様を調査しました。個人的にはwithin_apply関数の仕様が特に興味深いと感じました。回路最適化やシミュレーションの実行はサーバ側で走りますが、同時に複数プロセスから呼ぶとエラーになるようでした。これはバリエーション検証をする上で不便だと感じました。差別化の観点では、生成された回路サイズがqiskitより非常に小さくなることは印象的だと感じました。 qmodのベンダが保有する特許情報 8 を見ると、回路最適化の他にも誤り訂正やデバッグに関するものがあるようです。これらの要素が今後どのようにツールの機能に反映されていくのか興味を覚えました。