新しい量子素因数分解アルゴリズムの解説 (Part 5: SQIF の数値実験の分析)

本記事は, Yan 等の中国の研究者が 2022 年 12 月に発表した "[Factoring integers with sublinear resources on a superconducting quantum processor](https://arxiv.org/abs/2212.12372)" というプレプリントを解説した一連の記事の Part 5 です. この論文は QAOA という量子計算を使った新しい素因数分解アルゴリズムを提案しており, その名称を SQIF (Sublinear Quantum Integer Factorization) と呼んでいます. 本記事では SQIF の実験結果を分析します.

[Part 1: 論文の概要](https://qiita.com/izutetsuya/items/f40854aeb63d623410b8)
[Part 2: 平方差法](https://qiita.com/izutetsuya/items/0179d2609d666d95b392) 
[Part 3: Schnorr の素因数分解法](https://qiita.com/izutetsuya/items/a82aa95f11cf0de6ee6b)
[Part 4: SQIF の詳細](https://qiita.com/izutetsuya/items/2ceeb8c8dd581799b4da) 
Part 5: SQIF の実験結果の分析 (この記事)
Part 6: まとめ

https://arxiv.org/abs/2212.12372

# SQIF の処理手順

Yan et al. のプレプリントで紹介されている $N=1961$ の素因数分解実験を分析します. SQIF の処理手順は以下の通りです.
(1) 格子の生成
(2) 拡張関係式の収集 (CVP の求解)
(3) 行列計算
(4) $N$ の因数分解

# (1) 格子の生成

$N=1961$ のビット長は $m=\lceil \log_2{N} \rceil=11$ なので, 格子の次元は $n = \lfloor \frac{m}{\log_2{m}} \rfloor = \lfloor \frac{11}{\log_2{11}} = \lfloor 3.74 \rceil = 3$ となります. 従って, 素数基底 (拡張関係式の左辺で使用する素数の集合) は $n=3$ 個の素数の集合 $\mathbb{P}=\\{2,3,5\\}$ となります. 

各基底ベクトルの最後の行の成分は ($N^c \log_2{p_i}$ の代わりに) $\lceil 10^c \log_2{p_i} \rfloor$ を使用し, ウェイトは $c=1.5$ とします. また関数 $f$ は (ランダム置換の代わりに) $f(i)=\lceil i/2 \rfloor$ とします. このとき基底ベクトルは
```math
  \mathbb{b}_1 = \begin{pmatrix} f(1) \\ 0 \\ 0 \\ \lceil 10^c \log_2{p_1} \rfloor \end{pmatrix}
  = \begin{pmatrix} 1 \\ 0 \\ 0 \\ 22 \end{pmatrix}, \\
  \mathbb{b}_2 = \begin{pmatrix} 0 \\ f(2) \\ 0 \\ \lceil 10^c \log_2{p_2} \rfloor \end{pmatrix}
  = \begin{pmatrix} 0 \\ 1 \\ 0 \\ 35 \end{pmatrix}, \\
  \mathbb{b}_3 = \begin{pmatrix} 0 \\ 0 \\ f(3) \\ \lceil 10^c \log_2{p_3} \rfloor \end{pmatrix}
  = \begin{pmatrix} 0 \\ 0 \\ 2 \\ 51 \end{pmatrix}
```
となり, 格子
```math
  B = (\mathbb{b}_1, \mathbb{b}_2, \mathbb{b}_3) =
    \begin{pmatrix}
      1 & 0 & 0 \\
      0 & 1 & 0 \\
      0 & 0 & 2 \\
      22 & 35 & 51
    \end{pmatrix}
```
が得られました.  

# (2) 拡張関係式の収集 (CVP の求解)
## (2-1) ターゲットベクトルの設定

CVP を考える際のターゲットベクトルは
```math
  \mathbb{t} = \begin{pmatrix} 0 \\ 0 \\ 0 \\ \lceil 10^c \log_2{N} \rfloor \end{pmatrix}
  = \begin{pmatrix} 0 \\ 0 \\ 0 \\ 240 \end{pmatrix}
``` 
となります.

## (2-2) LLL アルゴリズムの適用

格子 $B$ に対して LLL アルゴリズムを適用して縮約格子
```math
  D = {\rm LLL}(B) =
    \begin{pmatrix}
       1 & -4 & -3 \\
      -2 &  1 &  2 \\
       2 &  2 &  0 \\
       3 & -2 &  4
    \end{pmatrix}
    = (\mathbb{d}_1,\mathbb{d}_2,\mathbb{d}_3)
```
を得ます. (プレプリントでは LLL アルゴリズムのパラメータが $\delta=0.75$ となっていますが, blueqat の湊さんの追試では列ベクトルが入れ替わってしまうようなので, パラメータは $\delta=0.99$ を使用したのではないかと思います.)

## (2-3) Babai アルゴリズムの適用

格子 $D$ に Babai アルゴリズムを適用してターゲットベクトル $\mathbb{t}$ に対する格子 $B$ にでの最近ベクトル問題の近似解
```math
  \mathbb{b}_{\rm op} = \begin{pmatrix} 0 \\ 4 \\ 4 \\ 242 \end{pmatrix}
```
を得ます. 念のため検算しておくと, この近似解は
```math
  \mathbb{b}_{\rm op} = 0 \mathbb{b}_1 + 4 \mathbb{b}_2 + 2  \mathbb{b}_3
```
と表せるため, 確かに CVP の近似解となる格子点であることわかります.

## (2-4) ハミルトニアンの構成

最適化問題を解くためのハミルトニアンベクトルは
```math
  \begin{eqnarray}
    \mathbb{h}
      &=& \mathbb{t} - \sum x_i \mathbb{d}_i  - \mathbb{b}_{\rm op} \\
      &=& \begin{pmatrix} 0 \\ 0 \\ 0 \\ 240 \end{pmatrix}
         - x_1 \begin{pmatrix} 1 \\ -2 \\ 2 \\ 3 \end{pmatrix} 
         - x_2 \begin{pmatrix} -4 \\ 1 \\ 2 \\ -2 \end{pmatrix} 
         - x_3 \begin{pmatrix} -3 \\ 2 \\ 0 \\ 4 \end{pmatrix} 
         - \begin{pmatrix} 0 \\ 4 \\ 4 \\ 242 \end{pmatrix} \\
      &=&  \begin{pmatrix}
             - x_1 + 4 x_2 + 3 x_3 \\
             2 x_1 - x_2 - 2 x_3 - 4 \\
             -2 x_1 - 2 x_2 - 4 \\
             -3 x_1 + 2 x_2 - 4 x_3 - 2
           \end{pmatrix}
  \end{eqnarray}
```
となります. よって最適化するハミルトニアンは
```math
  \begin{eqnarray}
    H(x_1,x_2,x_3) 
      &=& ||\mathbb{h}||^2 \\
      &=& (- x_1 + 4 x_2 + 3 x_3)^2
         + (2 x_1 - x_2 - 2 x_3 - 4)^2
         + (-2 x_1 - 2 x_2 - 4)^2
         + (-3 x_1 + 2 x_2 - 4 x_3 - 2)^2 \\
      &=& 36 + 30 x_1 + 41 x_2 + 61 x_3 - 16 x_1 x_2 + 12 x_2 x_3 + 10 x_3 x_1
  \end{eqnarray}
```
となります. ただし $x_1, x_2, x_3 \in \\{0,1\\}$ です.

## (2-5) 最適化問題の求解

次にこのハミルトニアンに対して QAOA などのソルバーを適用し, 最適化問題の解を求めます. われわれが扱っているハミルトニアン $H(x_1,x_2,x_3)$ は $2^3=8$ 通りの値しかとりませんので, ここではソルバーを使わずに全ての値を書き出して以下の表を得ます. ここで R3 は Babai のアルゴリズムによる近似解に対応しているので, それよりも良い解は R1, R2 となります.

| 名称 | $H$ | $x_1$ | $x_2$ | $x_3$ | $u$ | $v$ | $\| u-vN \|$ |
|:----:|:---:|:-----:|:-----:|:-----:|:--------------|:---:|:---------|
| R1   | $33$ | $0$ | $0$ | $1$ | $1800 = 2^3 \times 3^2 \times 5^2$ | $1$ | $7 \times 23$ |
| R2 | $35$ | $1$ | $1$ | $0$ | $1944 = 2^3 \times 3^5 \times 5^0$ | $1$ | $17$ |
| R3 | $36$ | $0$ | $0$ | $0$ | $2025 = 2^0 \times 3^4 \times 5^2$ | $1$ | $2^6$ |
| R4 | $42$ | $1$ | $0$ | $0$ | $3645 = 2^0 \times 3^6 \times 5^1$ | $2$ | $277$ |
| R5 | $45$ | $0$ | $1$ | $0$ | $2160 = 2^4 \times 3^3 \times 5^2$ | $1$ | $199$ |
| R6 | $49$ | $1$ | $0$ | $1$ | $1620 = 2^7 \times 3^1 \times 5^1$ | $1$ | $11 \times 31$ |
| R7 | $54$ | $0$ | $1$ | $1$ | $1920 = 2^7 \times 3^3 \times 5^1$ | $1$ | $41$ |
| R8 | $54$ | $1$ | $1$ | $1$ | $1728 = 2^6 \times 3^3 \times 5^0$ | $1$ | $233$ |

## (2-6) 拡張関係式の収集

プレプリントの実験で収集した拡張関係式は以下の通りです. 残念ながら週数方法は不明です.

| 名称 | $u$ | $v$ | $\| u-vN \|$ |
|:----:|:--------------|:---:|:---------|
| sn01 | $1944 = 2^3 \times 3^5 \times 5^0$ | $1$ | $17$ |
| sn02 | $2000 = 2^4 \times 3^0 \times 5^3$ | $1$ | $3 \times 13$ |
| sn03 | $1920 = 2^7 \times 3^1 \times 5^1$ | $1$ | $41$ |
| sn04 | $2025 = 2^0 \times 3^4 \times 5^2$ | $1$ | $2^6$ |
| sn05 | $1875 = 2^0 \times 3^1 \times 5^4$ | $1$ | $2 \times 43$ |
| sn06 | $1800 = 2^3 \times 3^2 \times 5^2$ | $1$ | $7 \times 23$ |
| sn07 | $2250 = 2^1 \times 3^2 \times 5^3$ | $1$ | $17^2$ |
| sn08 | $1620 = 2^2 \times 3^4 \times 5^1$ | $1$ | $11 \times 31$ |
| sn09 | $1600 = 2^6 \times 3^0 \times 5^2$ | $1$ | $19^2$ |
| sn10 | $2500 = 2^2 \times 3^0 \times 5^4$ | $1$ | $7^2 \times 11$ |
| sn11 | $1350 = 2^1 \times 3^3 \times 5^2$ | $1$ | $13 \times 47$ |
| sn12 | $1296 = 2^4 \times 3^4 \times 5^0$ | $1$ | $5 \times 7 \times 19$ |
| sn13 | $1125 = 2^0 \times 3^2 \times 5^3$ | $1$ | $2^2 \times 11 \times 19$ |
| sn14 | $1000 = 2^3 \times 3^0 \times 5^3$ | $1$ | $31^2$ |
| sn15 | $972 = 2^2 \times 3^5 \times 5^0$ | $1$ | $23 \times 43$ |
| sn16 | $800 = 2^5 \times 3^0 \times 5^2$ | $1$ | $3^3 \times 43$ |
| sn17 | $11664 = 2^4 \times 3^6 \times 5^0$ | $5$ | $11 \times 13^2$ |
| sn18 | $3888 = 2^4 \times 3^5 \times 5^0$ | $1$ | $41 \times 47$ |
| sn19 | $6075 = 2^0 \times 3^5 \times 5^2$ | $1$ | $2 \times 11^2 \times 17$ |
| sn20 | $9375 = 2^0 \times 3^1 \times 5^5$ | $2$ | $7 \times 19 \times 41$ |

ここで, sn01=R2, sn03=R7, sn04=R3, sn06=R1, sn8=R6 となっていますが, 残りの拡張関係式は別の方法で収集されています.

## (3) 行列計算

(2-6) で得られた20個の拡張関係式のうち, 16個は線型独立となっています. ここで拡張関係式で使用している素数は 13 個なので, 掃き出し法によって複数個の関係式を求められそうです. 詳細な計算は省略します.

## (4) 素因数分解

プレプリントでは素因数分解に成功する 4 通りの計算例が示されています. 以下, それぞれの計算内容を確認していきます.

## (4-1) Eg. 1

sn04 (=R3) の拡張関係式は全てのべき指数が偶数になっているので, それ自体が求めたい関係式になっています. そこで $x=3^2 \times 5 = 45$, $y=2^3=8$ とおけば, $\gcd(x-y,N)=\gcd(45-8,1961)=37$ が得られ, 素因数分解に成功します. ただしこの拡張関係式は Babai のアルゴリズムによって求められる解そのものなので, Yan et al. の提案法は寄与していないことに注意が必要です.

## (4-2) Eg. 2

sn09 の拡張関係式は全てのべき指数が偶数になっているので, それ自体が求めたい関係式になっています. $x=2^3 \times 5 = 40$, $y=19$ とおくと, $\gcd(x-y,N)=\gcd(40-19,1961)=1$ となり, 素因数分解に失敗します. (実はこの場合は $40^2 \equiv -19^2 \pmod{N}$ となっているため, 望ましい合同式になっていません).

## (4-3) Eg. 3

sn10 と sn17 を組み合わせることで
```math
  2^6 \times 3^6 \times 5^4 \equiv 7^2 \times 11^2 \times 17^2 \pmod{1961}
```
という拡張関係式が得られるので, $x=2^3 \times 3^3 \times 5^2=5400$, $y=7 \times 11 \times 17=1001$ とおくと, $\gcd(x-y,N)=\gcd(5400-1001,1961)=37$ となり, 素因数分解に成功します. ただしこれらの拡張関係式はハミルトニアンからは求められないので, Yan et al. の提案法は寄与していないことに注意が必要です.

## (4-4) Eg. 4

sn05 と sn16 を組み合わせることで
```math
  2^4 \times 5^6 \equiv 3^2 \times 43^2 \pmod{1961}
```
という拡張関係式が得られるので, $x=2^2 \times 5^3 =500$, $y=3 \times 43=129$ とおくと, $\gcd(x-y,N)=\gcd(500-129,1961)=53$ となり, 素因数分解に成功します. ただしこれらの拡張関係式はハミルトニアンからは求められないので, Yan et al. の提案法は寄与していないことに注意が必要です.

# まとめ

SQIF の計算例として, プレプリントで紹介されている $N=1961$ に対する具体的な計算処理を紹介しました. プレプリントでは 4 種類の組み合わせから素因数分解が求められると説明されていますが, そこで利用する拡張関係式はハミルトニアンからは得られないものであるため, 少なくともこの例では Yan et al. の提案手法の効果は不明となっています.