はじめに
こんにちは。量子コンピューティングを趣味でやっております。onewanです。量子コンピューター Advent Calendar 2025 の14日目に参加させていただいております。
さて今回は、量子コンピューターを用いたアルゴリズムとして有名なShorのアルゴリズムの中身を眺めて、自分なりに噛み砕いていきたいと思います。
Shorのアルゴリズムは、因数分解を高速で解くアルゴリズムで、理想的で大規模な量子コンピューターが完成すれば、Shorのアルゴリズムを利用することで、RSA暗号を高速で解くことが出来るといいます。
特定の条件下ではあるものの、100万量子ビット未満で1週間以内にRSA-2048を解読できるという推定がされているようです。
さて、今回は詳細な解説などは、既出の参考書やホームページを見ていただくとして、Shorのアルゴリズムの中身をざっくりと眺めて雰囲気を理解することが目的です。
Shorのアルゴリズムの実装
$15$の因数分解を行うShorのアルゴリズムについて、IBMのホームページでqiskitのコードを交えて解説されていますので、ご覧になってみてください。
概要
Shorのアルゴリズムは、ある合成数$N$を因数分解する手法です。
手順としては、まず$N$と互いに素な正整数$a(<N)$をとり、量子位相推定法を用いて$a$の位数を求め、位数を利用して因数分解を行います。
ここで出てきた用語について、それぞれ解説していきます。
- 位数
- 量子位相推定法
- 量子位相推定法を用いて位数を求める
- 位数を用いて因数分解を行う方法
位数
互いに素な$a,N \in \mathbb{Z},a<N$に対して$a^r \equiv 1 \pmod{N}$を満たす最小の正整数$r$のことを位数といいます。
$N$が素数であった場合にはフェルマーの小定理より、$a^{N-1} \equiv 1 \pmod{N}$となります。
$N$が素数とは限らない場合に拡張した定理として、オイラーの定理があり、$a^{\phi(N)} \equiv 1 \pmod{N}$が成立します。ここで、$\phi(N)$はオイラーのトーシェント関数と呼ばれるもので、$1$から$N-1$までで$N$と互いに素である数の個数を表します。
ただ、フェルマーの小定理（またはオイラーの定理）は、少なくとも$N-1$（または$\phi(N)$）乗すれば、$N$と互いに素な全ての$a$で$1$と同値になることを保証しているものであって、$a$を固定したときに、$N-1$乗して初めて$1$と同値になること（最小性）を保証するものでは無いことに注意して下さい。
実際、$a=10,N=37$のとき、$10^{r} \equiv 1 \pmod{37}$を満たす$r$としては、フェルマーの小定理から$10^{36} \equiv 1 \pmod{37}$が成立しますが、これは少なくとも$36$乗すれば$1$と同値になることを保証しているだけで、実際に最も小さい$r$は、その約数となります。実際、$10^{3} \equiv 1 \pmod{37}$となりますので、位数は$3$となります。
この性質は、トーシェント関数$\phi(N)$も同様で、オイラーの定理からは、位数がトーシェント関数$\phi(N)$の約数であることまでしか分かりません。
量子位相推定法
量子位相推定法は、多くの量子アルゴリズムで利用されるサブルーチンのようなもので、ユニタリー$U$と固有ベクトルから、対応する固有値を補助量子ビット数$m$に対して$2^m$精度で求める事が出来る手法です。
量子位相推定法を用いて位数を求める
以下の手順で進めます。
- 位数$r$を求めるためのユニタリー$M_a$を定義
- Phase estimation circuitの初期状態を準備するために、$M_a$の固有ベクトル$\ket{\psi}$を見つける
- $\frac{j}{r}$の近似値$\frac{y}{2^m}$を求める
$M_a \ket{x} = \ket{ax \pmod{N}}$となるように乗算演算子$M_a$を定義します。要は$a$倍して$\mod{N}$を取ったものです。ここで、$M_a$はユニタリー演算子になることが知られています。
任意の$j (\in \{ 0,\dots,r-1 \})$に対して、$\ket{\psi_j}=\frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} \omega_r^{-jk} \ket{a^k}$は$M_a$の固有ベクトルとなり、対応する固有値は$\omega^j_r = e^{2\pi i \frac{j}{r}}$となります。
実際、
$$M_a \ket{\psi_j}=\frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} \omega_r^{-jk} \ket{a^{k+1}}$$
$$\omega_r^j\ket{\psi_j}=\frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} \omega_r^{j(1-k)} \ket{a^k}$$
について、$\ket{a^x}$の係数を比較すると、上式は$k=x-1$の係数、下式は$k=x$の係数なので、それぞれ、（上式の$a^x$の係数）＝$\frac{1}{\sqrt{r}}\omega^{-j(x-1)}$で、（下式の$a^x$の係数）＝$\frac{1}{\sqrt{r}}\omega^{j(1-x)}$となり、$x=\{1,2,\dots,r-1 \}$で一致し、また$a^r \equiv 1=a^0 \pmod{N}$より$x=0$と$x=r$のときは一致することより、上式と下式は等しい事がわかります。
量子回路に入力するための固有ベクトルを、明に作ることは出来ませんが、$\ket{1}=\frac{1}{\sqrt{r}}\sum_{k=0}^{r-1} \ket{\psi_k}$なので、$\ket{1}$を入力してやると、$k\in\{0,\dots,r-1\}$を一様ランダムに選び、$\ket{\psi_k}$を入力した場合と同じ測定結果が得られます。
したがって、$k/r$ に対する近似値$\frac{y}{2^m}$の値が得られます。要は、これを繰り返せば、$0,\frac{1}{r},\frac{2}{r},\dots,\frac{r-1}{r}$が近似値として求まるので、そこから$\frac{1}{r}$の値が計算できて、$r$が求まるということです。
位数をもとに因数分解を行う
位数とは、互いに素な$a,N \in \mathbb{Z},a<N$に対して$a^r \equiv 1 \pmod{N}$を満たす最小の正整数$r$のことでした。
よって、ある正整数$k$があって、$a^r - 1= k N$となります。もし、$r$が偶数であれば、$(a^{\frac{r}{2}} - 1)(a^{\frac{r}{2}} + 1)= k N$となり、$N=GCD(N,(a^{\frac{r}{2}} - 1)) \times GCD(N,(a^{\frac{r}{2}} + 1))$となります。また、$r$が奇数のときは、新たな$a$を取り直して、位数$r$が偶数になるまで繰り返します。
さらに、$(a^{\frac{r}{2}} - 1)$ と$(a^{\frac{r}{2}} + 1)$にそれぞれ$N$の因数が１つずつ含まれていることがあります。そうすると、$GCD(N,(a^{\frac{r}{2}} - 1))$、もしくは$ GCD(N,(a^{\frac{r}{2}} + 1))$を求めることによって、$N$の因数を見つける事ができます。残念ながら$GCD$が$1$や$N$だった場合は、また$a$を取り直して繰り返します。
おわりに
そういうわけで、今回はShorのアルゴリズムについて、ざっくりと眺めて行きました。
余談ですが、量子位相推定法を用いずに、手計算で$r$を求めるということを某所で試みたのですが、たとえば$N=105$のとき、$105$と互いに素な$a$に対する位数を求めるのにオイラーの定理に出てくるトーシェントを求めようとしたのですが、$105$を素因数分解してしまえば簡単に求まるものの、それだと本末転倒なので、$105$と互いに素な数を数え上げる必要があり、おそらくユークリッドの互除法を使うことになると思いますが、$1$から$104$の数に対して手計算でそれを行うのがどれだけ大変かを考えると、位数を求める量子アルゴリズムのありがたさが分かりました。
あと、今回は触れなかったのですが、Qiskitのページで実装されているShorのアルゴリズムでは、途中でmodの計算を行うために、モジュラー指数演算子を構成しているのですが、具体的な例（N=15とa=2）でadhocなやり方で作成されているように見えたので、どうやって一般化するのか、またその構成にどのくらいの計算量が掛かるのかは、さらに調べてみたいと思っています。
それではみなさま、今年も残すところわずかとなりましたが、良い量子コンピューティングライフをお送りください。
その他参考サイト