最後に測定を行う。補助量子ビットを計測して $\ket{1}$ が観測されたケースに注目する（$\frac{1}{\lambda_j^\prime}$ を使いたいのでこちらに注目する）。この時（$\ket{0} \otimes \ket{0}^{\otimes n} \otimes \ket{b}$ に $QPE^{-1} U_f^\prime QPE$ を適用した結果の式）は、
$$
\begin{align*}
&\ \sum_{j=0}^{n-1} \left(\sqrt{1 - \left(\frac{C}{\lambda_j^\prime}\right)^2} \ket{0} + \frac{C}{\lambda_j^\prime} \ket{1} \right) \otimes \ket{0}^{\otimes n} \otimes \beta_j \ket{a_j} \\
=&\ C \sum_{j=0}^{n-1} \ket{1} \otimes \ket{0}^{\otimes n} \otimes \frac{\beta_j }{\lambda_j^\prime} \ket{a_j}
\tag{1}
\end{align*}
$$
となる。$C$ は正規化定数だったので、
$$
\begin{align*}
C = \sqrt{\sum_{j=0}^{n-1} \left|\frac{\beta_j }{\lambda_j^\prime}\right|^2}
\end{align*}
$$
である。
Qiskit (38) ― HHL アルゴリズム - らんだむな記憶 より
$$
\begin{align*}
\ket{x} = A^{-1} \ket{b} = \sum_{j=0}^{n-1} \frac{\beta_j }{\lambda_j} \ket{a_j}
\end{align*}
$$
を思い出すと、
$$
\begin{align*}
\sum_{j=0}^{n-1} \frac{\beta_j}{\lambda_j^\prime} \ket{a_j} = \frac{2 \pi} {t} \sum_{j=0}^{n-1} \frac{\beta_j}{\lambda_j} \ket{a_j} = \frac{2 \pi}{t} \ket{x}
\end{align*}
$$
であることが分かる。このことから (1) 式は
$$
\begin{align*}
\frac{2 \pi C}{t} \ket{1} \otimes \ket{0}^{\otimes n} \otimes \ket{x}
\end{align*}
$$
となる。よって、第 3 レジスタを測定して得た “バイナリ表示” から整数 $x^\prime$ を求め（例えば $\ket{101}$ が得られたなら $x^\prime = 2^2 \times 1 + 2^0 = 5$ とする）、$\frac{t}{2 \pi C}$ をかけ、最後に一番最初の正規化を解除するために $|b|$ をかけることで、$Ax = b$ の解 $x$ が求まる。
ここまでのまとめ
全体を俯瞰すると、要するに $U_f^\prime$ の操作でやりたかったのは、$CR_y(\theta)$ の適用で、位相 $\theta$ を工夫して補助レジスタに $A$ の固有値の逆数を取り出すことであったことが分かる。こうすることで、他の部分の情報と組み合わせてエルミート行列 $A$ のスペクトル分解から $A^{-1}$ が現れてくることになる。回転ゲートの位相 $\theta$ をどうやって決めているのかいまいちよく分からないのだが、今回はちょっと気にしないでおく。
ここまでご都合主義で話を単純化してもこの難易度なので、HHL アルゴリズムは相当難しい内容だなと思う。暗黙的には以下のようなことが実行されていることになる:
- $\ket{b}$ の固有ベクトルによる展開を考える（具体的な固有ベクトルは求めない）
- 量子位相推定で固有ベクトルに対する固有値（の位相）を求める
- ユニタリ操作で固有値の逆数を求める
- $\ket{b}$ の固有ベクトルによる展開と $A$ の固有値の逆数からスペクトル分解を通じて $A^{-1}$ を構成し、展開係数と $A^{-1}$ から $A^{-1} \ket{b}$ として $\ket{x}$ を求める
- 正規化を解除して連立方程式の解 $x$ を求める
—————
書籍の付属コードを見ると、若干変数名とかが分かりにくいが、Solving Linear Systems of Equations using HHL のコードを見ながら同ページ内の記号類を見ると、何を意味しているかはより明確になると思う。qrb が第 3 レジスタで $\ket{b}$ とかを入れている箇所で、qrl が第 2 レジスタで $\lambda_j^\prime$ とかを入れている箇所で、qra が補助ビットである。a real quantum device 向けの textbook コードに近いが、値として最適値を使うようにしてより簡単にしたような内容になっていると思われる。
import numpy as np A = np.array([[1, -1/3], [-1/3, 1]]) b = np.array([1, 0]) np.linalg.solve(A, b)
array([1.125, 0.375])
なのでこれに近い値が求まれば良い。補助ビット qra の値を古典レジスタ 0 で測定するように実装されているので、$\ket{...1}$ の観測を見ると良い。このケースが $\frac{1}{3} \ket{x} = p \ket{0} + q \ket{1}$ の重ね合わせであるということなので、ベクトル解を格納するレジスタ qrb が古典レジスタ 3 で測定されることから、$\ket{0..1} = \ket{0001}$ の観測確率が $\ket{0}$ の確率振幅 $p$ に、$\ket{1..1} = \ket{1001}$ の観測確率が $\ket{1}$ の確率振幅 $q$ に対応する。試しに実行したケースでは
import math 3*math.sqrt(0.136), 3*math.sqrt(0.017)
で
(1.1063453348751464, 0.39115214431215894)
なのでまぁまぁ近いと言えるだろう・・・。多分論文 [0811.3171] Quantum algorithm for solving linear systems of equations をちゃんと読まないと理解できないと思うので*1、雰囲気だけ感じたということにして p.159 まで完了とする。
qiskit-terra/qiskit/algorithms/linear_solvers/hhl.py#L466-L473 の感じだと何かしらうまくできるのだろう・・・。qiskit-terra/qiskit/algorithms/linear_solvers/hhl.py#L441-L443 の辺りにヒントがありそう・・・。textbook の図の中に Repeat untill success という文章が見えるので、幾らかイテレーションを回すのか・・・？
*1:該当箇所周辺を読んでみたけど、Adding an ancilla qubit and rotating conditioned on $\tilde{\lambda}_k$ yields で流されるので、やはり Qiskit とかで一般的にどう実装するのかよく分からなかった・・・
