쇼어의 알고리즘

사용량 추정치: Eagle r3 프로세서에서 3초 (참고: 이는 추정치일 뿐입니다. 런타임은 다를 수 있습니다.)

1994년 피터 쇼어가 개발한 쇼 알고리즘은 다항식 시간에서 정수를 인수분해하는 획기적인 양자 알고리즘입니다. 이 알고리즘의 중요성은 큰 정수를 알려진 기존 알고리즘보다 기하급수적으로 빠르게 인수분해할 수 있다는 점으로, 큰 수를 인수분해하는 데 어려움을 겪는 RSA와 같이 널리 사용되는 암호화 시스템의 보안을 위협할 수 있다는 것입니다. 쇼의 알고리즘은 충분히 강력한 양자 컴퓨터에서 이 문제를 효율적으로 해결함으로써 암호화, 사이버 보안, 계산 수학과 같은 분야에 혁명을 일으켜 양자 계산의 혁신적 힘을 강조할 수 있습니다.

이 튜토리얼에서는 양자 컴퓨터에서 15를 인수분해하여 쇼의 알고리즘을 시연하는 데 중점을 둡니다.

먼저 순서 찾기 문제를 정의하고 양자 위상 추정 프로토콜에서 해당 회로를 구성합니다. 다음으로, 트랜스파일할 수 있는 최단 심도의 회로를 사용하여 실제 하드웨어에서 순서 찾기 회로를 실행합니다. 마지막 섹션에서는 순서 찾기 문제를 정수 인수분해에 연결하여 쇼의 알고리즘을 완성합니다.

일반적인 구현과 15, 21과 같은 특정 정수를 인수분해하는 데 초점을 맞춰 실제 하드웨어에서 쇼 알고리즘의 다른 데모에 대해 논의하는 것으로 튜토리얼을 마무리합니다.

참고: 이 튜토리얼에서는 쇼의 알고리즘과 관련된 회로의 구현과 데모에 더 중점을 둡니다. 해당 자료에 대한 심층적인 교육 자료를 보려면 John Watrous 박사의 양자 알고리즘 기초 과정과 참고문헌 섹션의 논문을 참조하세요.

요구사항

이 튜토리얼을 시작하기 전에 다음이 설치되어 있는지 확인하세요:

Qiskit SDK v2.0 또는 이후 버전, 시각화 지원 기능 포함
Qiskit Runtime v0.40 이상 (pip install qiskit-ibm-runtime)

설정

import numpy as np
import pandas as pd
from fractions import Fraction
from math import floor, gcd, log

from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister
from qiskit.circuit.library import QFT, UnitaryGate
from qiskit.transpiler import CouplingMap, generate_preset_pass_manager
from qiskit.visualization import plot_histogram

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import SamplerV2 as Sampler

1단계: 고전적 입력을 양자 문제에 매핑하기

배경

쇼의 정수 인수분해 알고리즘은 순서 찾기 문제로 알려진 중간 문제를 활용합니다. 이 섹션에서는 양자 위상 추정을 사용하여 순서 찾기 문제를 해결하는 방법을 설명합니다.

위상 추정 문제

위상 추정 문제에서는 $n$ 큐비트의 양자 상태 $\ket{\psi}$ 와 $n$ 큐비트에 작용하는 단일 양자 회로가 주어집니다. $\ket{\psi}$ 은 회로의 동작을 설명하는 행렬 $U$ 의 고유 벡터이며, 우리의 목표는 $\ket{\psi}$ 에 해당하는 고유값 $\lambda = e^{2 \pi i \theta}$ 을 계산하거나 근사값을 구하는 것입니다. 즉, 회로는 다음을 만족하는 $\theta \in [0, 1)$ 수에 대한 근사값을 출력해야 합니다 $U \ket{\psi}= e^{2 \pi i \theta} \ket{\psi}.$ 위상 추정 회로의 목표는 $\theta$ 을 $m$ 비트로 근사화하는 것입니다. 수학적으로 말하자면, $\theta \approx y / 2^m$ 이 $y \in {0, 1, 2, \dots, 2^{m-1}}$ 이 되도록 $y$ 을 구하고 싶습니다. 다음 이미지는 $m$ 큐비트에서 측정하여 $m$ 비트에서 $y$ 을 추정하는 양자 회로를 보여줍니다.

위 회로에서 상위 $m$ 큐비트는 $\ket{0^m}$ 상태에서 초기화되고, 하위 $n$ 큐비트는 $\ket{\psi}$ 상태에서 초기화됩니다. 이는 $U$ 의 고유벡터가 됩니다. 위상 추정 회로의 첫 번째 구성 요소는 해당 제어 큐비트로 위상 킥백을 수행하는 제어 단위 연산입니다. 이러한 제어 부호는 최하위 비트부터 최상위 비트까지 제어 큐비트의 위치에 따라 지수화됩니다. $\ket{\psi}$ 는 $U$ 의 고유 벡터이므로 아래쪽 $n$ 큐비트의 상태는 이 연산에 영향을 받지 않지만, 고유값의 위상 정보는 위쪽 $m$ 큐비트로 전파됩니다.

제어 유니티를 통한 위상 킥백 연산 후, 상위 $m$ 큐비트의 모든 가능한 상태는 유니티 $U$ 의 각 고유 벡터 $\ket{\psi}$ 에 대해 서로 직교하는 것으로 나타났습니다. 따라서 이러한 상태는 완벽하게 구분할 수 있으며, 이러한 상태가 형성하는 기저를 다시 계산 기저로 회전하여 측정할 수 있습니다. 수학적 분석에 따르면 이 회전 행렬은 $2^m$ -차원 힐베르트 공간에서 역 양자 푸리에 변환(QFT)에 해당합니다. 모듈식 지수화 연산자의 주기적 구조가 양자 상태로 인코딩되고, QFT는 이 주기성을 주파수 영역에서 측정 가능한 피크로 변환한다는 것이 직관적입니다.

쇼의 알고리즘에 QFT 회로가 사용되는 이유에 대해 더 깊이 이해하려면 양자 알고리즘의 기초 과정을 참조하시기 바랍니다.

이제 위상 추정 회로를 순서 찾기에 사용할 준비가 되었습니다.

주문 찾기 문제

순서 찾기 문제를 정의하기 위해 몇 가지 수 이론 개념으로 시작합니다. 먼저, 주어진 양의 정수 $N$ 에 대해 $\mathbb{Z}_N$ 집합을 다음과 같이 정의합니다 $\mathbb{Z}_N = \{0, 1, 2, \dots, N-1\}.$ $\mathbb{Z}_N$ 의 모든 산술 연산은 모듈로 $N$ 로 수행됩니다. 특히 $N$ 과 공범인 모든 원소 $a \in \mathbb{Z}_n$ 는 특수하며 $\mathbb{Z}^*_N$ 를 다음과 같이 구성합니다 $\mathbb{Z}^*_N = \{ a \in \mathbb{Z}_N : \mathrm{gcd}(a, N)=1 \}.$ 원소 $a \in \mathbb{Z}^*_N$ 의 경우 가장 작은 양의 정수 $r$ 는 $a^r \equiv 1 \; (\mathrm{mod} \; N)$ 가 $a$ 모듈로 $N$ 의 순서로 정의됩니다. 나중에 살펴 보겠지만 $a \in \mathbb{Z}^*_N$ 의 순서를 찾으면 $N$ 을 인수 분해 할 수 있습니다.

위상 추정 회로에서 순서 찾기 회로를 구성하려면 두 가지 고려 사항이 필요합니다. 먼저, $r$ 차수를 구할 수 있는 유니타리 $U$ 를 정의하고, 둘째, 위상 추정 회로의 초기 상태를 준비하기 위해 $U$ 의 고유 벡터 $\ket{\psi}$ 를 정의해야 합니다.

순서 찾기 문제를 위상 추정과 연결하기 위해, 고전적 상태가 $\mathbb{Z}_N$ 에 해당하는 시스템에서 정의된 연산에 고정 요소 $a \in \mathbb{Z}^*_N$ 를 곱하는 연산을 고려합니다. 특히, 우리는 이 곱셈 연산자 $M_a$ 를 정의하여 $x \in \mathbb{Z}_N$ 각각에 대해 $M_a \ket{x} = \ket{ax \; (\mathrm{mod} \; N)}$ 이 되도록 합니다. 방정식의 오른쪽에 있는 케트 내부에서 곱셈 모듈로 $N$ 를 취하고 있음을 암시합니다. 수학적 분석에 따르면 $M_a$ 은 단일 연산자입니다. 또한 $M_a$ 에는 $a$ 의 순서 $r$ 를 위상 추정 문제에 연결할 수 있는 고유벡터와 고유값 쌍이 있다는 것을 알 수 있습니다. 구체적으로, $j \in \{0, \dots, r-1\}$ 의 선택에 대해 $\ket{\psi_j} = \frac{1}{\sqrt{r}} \sum^{r-1}_{k=0} \omega^{-jk}_{r} \ket{a^k}$ 은 $M_a$ 의 고유 벡터이며 해당 고유값은 $\omega^{j}_{r}$ 입니다 $\omega^{j}_{r} = e^{2 \pi i \frac{j}{r}}.$

관찰을 통해 편리한 고유벡터/고윳값 쌍은 $\ket{\psi_1}$ 와 $\omega^{1}_{r} = e^{2 \pi i \frac{1}{r}}$ 상태임을 알 수 있습니다. 따라서 고유 벡터 $\ket{\psi_1}$ 를 구할 수 있다면 양자 회로로 $\theta=1/r$ 위상을 추정할 수 있고, 따라서 $r$ 차수의 추정치를 얻을 수 있습니다. 그러나 그렇게 하는 것은 쉽지 않으므로 대안을 고려해야 합니다.

계산 상태 $\ket{1}$ 를 초기 상태로 준비하면 회로가 어떤 결과를 가져올지 생각해 보겠습니다. 이것은 $M_a$ 의 고유 상태는 아니지만 위에서 설명한 고유 상태의 균일한 중첩입니다. 즉, 다음과 같은 관계가 성립합니다. $\ket{1} = \frac{1}{\sqrt{r}} \sum^{r-1}_{k=0} \ket{\psi_k}$ 위 방정식의 의미는 초기 상태를 $\ket{1}$ 로 설정하면 위상 추정 회로에서 $k \in \{ 0, \dots, r-1\}$ 을 무작위로 선택하고 $\ket{\psi_k}$ 을 고유 벡터로 사용한 경우와 정확히 동일한 측정 결과를 얻을 수 있다는 것입니다. 즉, 상위 $m$ 큐비트를 측정하면 $k / r$ 값에 대한 근사치 $y / 2^m$ 가 산출되며, 여기서 $k \in \{ 0, \dots, r-1\}$ 은 무작위로 균일하게 선택됩니다. 이를 통해 저희의 목표였던 여러 번의 독립 실행 끝에 높은 수준의 자신감을 가지고 $r$ 을 학습할 수 있었습니다.

모듈러 지수 연산자

지금까지 양자 회로에서 $U = M_a$ 및 $\ket{\psi} = \ket{1}$ 을 정의하여 위상 추정 문제를 순서 찾기 문제와 연결했습니다. 따라서 마지막으로 남은 요소는 $M_a$ 의 모듈 지수를 $k = 1, 2, 4, \dots, 2^{m-1}$ 에 대해 $M_a^k$ 로 정의하는 효율적인 방법을 찾는 것입니다. 이 계산을 수행하기 위해 선택한 모든 거듭제곱 $k$ 에 대해 $M_a^k$ 에 대한 회로를 $M_a$ 에 대한 회로의 $k$ 배를 반복하는 대신 $b = a^k \; \mathrm{mod} \; N$ 를 계산한 다음 $M_b$ 에 대한 회로를 사용하면 된다는 것을 알 수 있습니다. 2의 거듭제곱만 필요하므로 반복 제곱을 사용하면 고전적으로 효율적으로 이 작업을 수행할 수 있습니다.

2단계: 양자 하드웨어 실행을 위한 문제 최적화

$N = 15$ 및 $a=2$ 의 구체적인 예시

여기서 잠시 멈추고 구체적인 예를 살펴보고 $N=15$ 에 대한 순서 찾기 회로를 구성해 보겠습니다. $N=15$ 에 대해 가능한 $a \in \mathbb{Z}_N^*$ 는 $a \in \{2, 4, 7, 8, 11, 13, 14 \}$ 입니다. 이 예제에서는 $a=2$ 를 선택합니다. $M_2$ 연산자와 모듈식 지수 연산자 $M_2^k$ 를 구성하겠습니다.

계산 기준 상태에서 $M_2$ 의 동작은 다음과 같습니다. $M_2 \ket{0} = \ket{0} \quad M_2 \ket{5} = \ket{10} \quad M_2 \ket{10} = \ket{5}$ $M_2 \ket{1} = \ket{2} \quad M_2 \ket{6} = \ket{12} \quad M_2 \ket{11} = \ket{7}$ $M_2 \ket{2} = \ket{4} \quad M_2 \ket{7} = \ket{14} \quad M_2 \ket{12} = \ket{9}$ $M_2 \ket{3} = \ket{6} \quad M_2 \ket{8} = \ket{1} \quad M_2 \ket{13} = \ket{11}$ $M_2 \ket{4} = \ket{8} \quad M_2 \ket{9} = \ket{3} \quad M_2 \ket{14} = \ket{13}$ 관찰을 통해 기저 상태가 섞여 있으므로 순열 행렬이 있다는 것을 알 수 있습니다. 스왑 게이트가 있는 4개의 큐비트에서 이 연산을 구성할 수 있습니다. 아래에서는 $M_2$ 및 제어- $M_2$ 작업을 구성합니다.

def M2mod15():
    """
    M2 (mod 15)
    """
    b = 2
    U = QuantumCircuit(4)

    U.swap(2, 3)
    U.swap(1, 2)
    U.swap(0, 1)

    U = U.to_gate()
    U.name = f"M_{b}"

    return U

# Get the M2 operator
M2 = M2mod15()

# Add it to a circuit and plot
circ = QuantumCircuit(4)
circ.compose(M2, inplace=True)
circ.decompose(reps=2).draw(output="mpl", fold=-1)

Output:

def controlled_M2mod15():
    """
    Controlled M2 (mod 15)
    """
    b = 2
    U = QuantumCircuit(4)

    U.swap(2, 3)
    U.swap(1, 2)
    U.swap(0, 1)

    U = U.to_gate()
    U.name = f"M_{b}"
    c_U = U.control()

    return c_U

# Get the controlled-M2 operator
controlled_M2 = controlled_M2mod15()

# Add it to a circuit and plot
circ = QuantumCircuit(5)
circ.compose(controlled_M2, inplace=True)
circ.decompose(reps=1).draw(output="mpl", fold=-1)

Output:

두 개 이상의 큐비트에 작용하는 게이트는 두 개의 큐비트 게이트로 더 세분화됩니다.

circ.decompose(reps=2).draw(output="mpl", fold=-1)

Output:

이제 모듈식 지수 연산자를 구성해야 합니다. 위상 추정에서 충분한 정밀도를 얻기 위해 추정 측정에 8개의 큐비트를 사용합니다. 따라서 각 $k = 0, 1, \dots, 7$ 에 대해 $b = a^{2^k} \; (\mathrm{mod} \; N)$ 을 사용하여 $M_b$ 을 구성해야 합니다.

def a2kmodN(a, k, N):
    """Compute a^{2^k} (mod N) by repeated squaring"""
    for _ in range(k):
        a = int(np.mod(a**2, N))
    return a

k_list = range(8)
b_list = [a2kmodN(2, k, 15) for k in k_list]

print(b_list)

Output:

[2, 4, 1, 1, 1, 1, 1, 1]

$b$ 값 목록에서 볼 수 있듯이 이전에 구축한 $M_2$ 외에도 $M_4$ 및 $M_1$ 도 구축해야 합니다. $M_1$ 은 계산 기반 상태에서 사소한 역할을 하므로 단순히 ID 연산자에 불과합니다.

$M_4$ 는 다음과 같이 계산 기반 상태에서 작동합니다. $M_4 \ket{0} = \ket{0} \quad M_4 \ket{5} = \ket{5} \quad M_4 \ket{10} = \ket{10}$ $M_4 \ket{1} = \ket{4} \quad M_4 \ket{6} = \ket{9} \quad M_4 \ket{11} = \ket{14}$ $M_4 \ket{2} = \ket{8} \quad M_4 \ket{7} = \ket{13} \quad M_4 \ket{12} = \ket{3}$ $M_4 \ket{3} = \ket{12} \quad M_4 \ket{8} = \ket{2} \quad M_4 \ket{13} = \ket{7}$ $M_4 \ket{4} = \ket{1} \quad M_4 \ket{9} = \ket{6} \quad M_4 \ket{14} = \ket{11}$

따라서 이 순열은 다음과 같은 스왑 연산으로 구성할 수 있습니다.

def M4mod15():
    """
    M4 (mod 15)
    """
    b = 4
    U = QuantumCircuit(4)

    U.swap(1, 3)
    U.swap(0, 2)

    U = U.to_gate()
    U.name = f"M_{b}"

    return U

# Get the M4 operator
M4 = M4mod15()

# Add it to a circuit and plot
circ = QuantumCircuit(4)
circ.compose(M4, inplace=True)
circ.decompose(reps=2).draw(output="mpl", fold=-1)

Output:

def controlled_M4mod15():
    """
    Controlled M4 (mod 15)
    """
    b = 4
    U = QuantumCircuit(4)

    U.swap(1, 3)
    U.swap(0, 2)

    U = U.to_gate()
    U.name = f"M_{b}"
    c_U = U.control()

    return c_U

# Get the controlled-M4 operator
controlled_M4 = controlled_M4mod15()

# Add it to a circuit and plot
circ = QuantumCircuit(5)
circ.compose(controlled_M4, inplace=True)
circ.decompose(reps=1).draw(output="mpl", fold=-1)

Output:

두 개 이상의 큐비트에 작용하는 게이트는 두 개의 큐비트 게이트로 더 세분화됩니다.

circ.decompose(reps=2).draw(output="mpl", fold=-1)

Output:

주어진 $b \in \mathbb{Z}^*_N$ 에 대한 $M_b$ 연산자는 순열 연산이라는 것을 확인했습니다. 여기서 다루는 순열 문제의 크기가 상대적으로 작기 때문에 $N=15$ 에는 4개의 큐비트만 필요하기 때문에 검사를 통해 SWAP 게이트로 이러한 연산을 직접 합성할 수 있었습니다. 일반적으로 이는 확장 가능한 접근 방식이 아닐 수 있습니다. 대신, 순열 행렬을 명시적으로 구성하고 키스킷의 UnitaryGate 클래스와 트랜스필레이션 메서드를 사용하여 이 순열 행렬을 합성해야 할 수도 있습니다. 하지만 이렇게 하면 회로가 훨씬 더 깊어질 수 있습니다. 예제는 다음과 같습니다.

def mod_mult_gate(b, N):
    """
    Modular multiplication gate from permutation matrix.
    """
    if gcd(b, N) > 1:
        print(f"Error: gcd({b},{N}) > 1")
    else:
        n = floor(log(N - 1, 2)) + 1
        U = np.full((2**n, 2**n), 0)
        for x in range(N):
            U[b * x % N][x] = 1
        for x in range(N, 2**n):
            U[x][x] = 1
        G = UnitaryGate(U)
        G.name = f"M_{b}"
        return G

# Let's build M2 using the permutation matrix definition
M2_other = mod_mult_gate(2, 15)

# Add it to a circuit
circ = QuantumCircuit(4)
circ.compose(M2_other, inplace=True)
circ = circ.decompose()

# Transpile the circuit and get the depth
coupling_map = CouplingMap.from_line(4)
pm = generate_preset_pass_manager(coupling_map=coupling_map)
transpiled_circ = pm.run(circ)

print(f"qubits: {circ.num_qubits}")
print(
    f"2q-depth: {transpiled_circ.depth(lambda x: x.operation.num_qubits==2)}"
)
print(f"2q-size: {transpiled_circ.size(lambda x: x.operation.num_qubits==2)}")
print(f"Operator counts: {transpiled_circ.count_ops()}")
transpiled_circ.decompose().draw(
    output="mpl", fold=-1, style="clifford", idle_wires=False
)

Output:

qubits: 4
2q-depth: 94
2q-size: 96
Operator counts: OrderedDict({'cx': 45, 'swap': 32, 'u': 24, 'u1': 7, 'u3': 4, 'unitary': 3, 'circuit-335': 1, 'circuit-338': 1, 'circuit-341': 1, 'circuit-344': 1, 'circuit-347': 1, 'circuit-350': 1, 'circuit-353': 1, 'circuit-356': 1, 'circuit-359': 1, 'circuit-362': 1, 'circuit-365': 1, 'circuit-368': 1, 'circuit-371': 1, 'circuit-374': 1, 'circuit-377': 1, 'circuit-380': 1})

이 개수를 수동으로 구현한 $M_2$ 게이트의 컴파일된 회로 깊이와 비교해 보겠습니다.

# Get the M2 operator from our manual construction
M2 = M2mod15()

# Add it to a circuit
circ = QuantumCircuit(4)
circ.compose(M2, inplace=True)
circ = circ.decompose(reps=3)

# Transpile the circuit and get the depth
coupling_map = CouplingMap.from_line(4)
pm = generate_preset_pass_manager(coupling_map=coupling_map)
transpiled_circ = pm.run(circ)

print(f"qubits: {circ.num_qubits}")
print(
    f"2q-depth: {transpiled_circ.depth(lambda x: x.operation.num_qubits==2)}"
)
print(f"2q-size: {transpiled_circ.size(lambda x: x.operation.num_qubits==2)}")
print(f"Operator counts: {transpiled_circ.count_ops()}")
transpiled_circ.draw(
    output="mpl", fold=-1, style="clifford", idle_wires=False
)

Output:

qubits: 4
2q-depth: 9
2q-size: 9
Operator counts: OrderedDict({'cx': 9})

보시다시피 순열 행렬 접근 방식은 단일 $M_2$ 게이트에 대해서도 수동 구현에 비해 훨씬 더 깊은 회로를 생성했습니다. 따라서 기존 운영 방식( $M_b$ )을 계속 유지할 예정입니다.

이제 앞서 정의한 제어 모듈식 지수 연산자를 사용하여 전체 차수 찾기 회로를 구성할 준비가 되었습니다. 다음 코드에서는 각 큐비트에 Hadamard 게이트, 일련의 controlled-U1 (또는 위상에 따라 Z) 게이트, 스왑 게이트 레이어를 사용하는 Qiskit Circuit 라이브러리에서 QFT 회로를 가져옵니다.

# Order finding problem for N = 15 with a = 2
N = 15
a = 2

# Number of qubits
num_target = floor(log(N - 1, 2)) + 1  # for modular exponentiation operators
num_control = 2 * num_target  # for enough precision of estimation

# List of M_b operators in order
k_list = range(num_control)
b_list = [a2kmodN(2, k, 15) for k in k_list]

# Initialize the circuit
control = QuantumRegister(num_control, name="C")
target = QuantumRegister(num_target, name="T")
output = ClassicalRegister(num_control, name="out")
circuit = QuantumCircuit(control, target, output)

# Initialize the target register to the state |1>
circuit.x(num_control)

# Add the Hadamard gates and controlled versions of the
# multiplication gates
for k, qubit in enumerate(control):
    circuit.h(k)
    b = b_list[k]
    if b == 2:
        circuit.compose(
            M2mod15().control(), qubits=[qubit] + list(target), inplace=True
        )
    elif b == 4:
        circuit.compose(
            M4mod15().control(), qubits=[qubit] + list(target), inplace=True
        )
    else:
        continue  # M1 is the identity operator

# Apply the inverse QFT to the control register
circuit.compose(QFT(num_control, inverse=True), qubits=control, inplace=True)

# Measure the control register
circuit.measure(control, output)

circuit.draw("mpl", fold=-1)

Output:

$M_1$ 은 동일성 연산자이므로 나머지 제어 큐비트에서 제어 모듈식 지수 연산을 생략했습니다.

이 튜토리얼의 뒷부분에서는 ibm_marrakesh 백엔드에서 이 회로를 실행할 것입니다. 이를 위해 이 특정 백엔드에 따라 회로를 트랜스파일하고 회로 깊이와 게이트 수를 보고합니다.

service = QiskitRuntimeService()
backend = service.backend("ibm_marrakesh")
pm = generate_preset_pass_manager(optimization_level=2, backend=backend)

transpiled_circuit = pm.run(circuit)

print(
    f"2q-depth: {transpiled_circuit.depth(lambda x: x.operation.num_qubits==2)}"
)
print(
    f"2q-size: {transpiled_circuit.size(lambda x: x.operation.num_qubits==2)}"
)
print(f"Operator counts: {transpiled_circuit.count_ops()}")
transpiled_circuit.draw(
    output="mpl", fold=-1, style="clifford", idle_wires=False
)

Output:

2q-depth: 187
2q-size: 260
Operator counts: OrderedDict({'sx': 521, 'rz': 354, 'cz': 260, 'measure': 8, 'x': 4})

3단계: `Qiskit primitives` 명령어로 실행합니다

먼저 이상적인 시뮬레이터에서 이 회로를 실행하면 이론적으로 얻을 수 있는 결과에 대해 논의합니다. 아래는 1024개의 샷을 사용한 위 회로의 시뮬레이션 결과입니다. 보시다시피, 제어 큐비트에 대한 4개의 비트 문자열에 대해 거의 균일한 분포를 얻습니다.

# Obtained from the simulator
counts = {"00000000": 264, "01000000": 268, "10000000": 249, "11000000": 243}

plot_histogram(counts)

Output:

제어 큐비트를 측정하여 $M_a$ 연산자의 8비트 위상 추정값을 구합니다. 이 이진 표현을 10진수로 변환하여 측정된 위상을 찾을 수 있습니다. 위의 히스토그램에서 볼 수 있듯이 네 개의 서로 다른 비트스트링이 측정되었으며, 각 비트스트링은 다음과 같이 위상 값에 해당합니다.

# Rows to be displayed in table
rows = []
# Corresponding phase of each bitstring
measured_phases = []

for output in counts:
    decimal = int(output, 2)  # Convert bitstring to decimal
    phase = decimal / (2**num_control)  # Find corresponding eigenvalue
    measured_phases.append(phase)
    # Add these values to the rows in our table:
    rows.append(
        [
            f"{output}(bin) = {decimal:>3}(dec)",
            f"{decimal}/{2 ** num_control} = {phase:.2f}",
        ]
    )

# Print the rows in a table
headers = ["Register Output", "Phase"]
df = pd.DataFrame(rows, columns=headers)
print(df)

Output:

            Register Output           Phase
0  00000000(bin) =   0(dec)    0/256 = 0.00
1  01000000(bin) =  64(dec)   64/256 = 0.25
2  10000000(bin) = 128(dec)  128/256 = 0.50
3  11000000(bin) = 192(dec)  192/256 = 0.75

측정된 모든 위상은 $\theta = k / r$ 에 해당하며 $k$ 은 $\{0, 1, \dots, r-1 \}$ 에서 무작위로 균일하게 샘플링됩니다. 따라서 연속 분수 알고리즘을 사용하여 $k$ 을 구하고 $r$ 을 구할 수 있습니다. Python 에는 이 기능이 내장되어 있습니다. 예를 들어 fractions 모듈을 사용하여 플로트를 Fraction 객체로 바꿀 수 있습니다:

Fraction(0.666)

Output:

Fraction(5998794703657501, 9007199254740992)

이렇게 하면 결과를 정확히 반환하는 분수(이 경우 0.6660000...)가 제공되므로 위와 같은 엉뚱한 결과가 나올 수 있습니다. .limit_denominator() 메서드를 사용하여 분모가 특정 값 이하인 플로트와 가장 유사한 분수를 구할 수 있습니다:

# Get fraction that most closely resembles 0.666
# with denominator < 15
Fraction(0.666).limit_denominator(15)

Output:

Fraction(2, 3)

훨씬 낫습니다. 차수(r)는 N보다 작아야 하므로 최대 분모를 15 로 설정합니다:

# Rows to be displayed in a table
rows = []

for phase in measured_phases:
    frac = Fraction(phase).limit_denominator(15)
    rows.append(
        [phase, f"{frac.numerator}/{frac.denominator}", frac.denominator]
    )

# Print the rows in a table
headers = ["Phase", "Fraction", "Guess for r"]
df = pd.DataFrame(rows, columns=headers)
print(df)

Output:

   Phase Fraction  Guess for r
0   0.00      0/1            1
1   0.25      1/4            4
2   0.50      1/2            2
3   0.75      3/4            4

측정된 고유값 중 두 개가 올바른 결과를 제공했음을 알 수 있습니다: $r=4$ 이며, Shor의 순서 찾기 알고리즘이 실패할 가능성이 있음을 알 수 있습니다. 이러한 나쁜 결과는 $k = 0$ 또는 $k$ 과 $r$ 이 공범이 아니고 $r$ 대신 $r$ 의 계수가 주어졌기 때문입니다. 이에 대한 가장 쉬운 해결책은 $r$ 에 대한 만족스러운 결과를 얻을 때까지 실험을 반복하는 것입니다.

지금까지 시뮬레이터에서 위상 추정 회로를 사용하여 $a=2$ 으로 $N=15$ 의 순서 찾기 문제를 구현했습니다. 쇼 알고리즘의 마지막 단계는 순서 찾기 문제를 정수 인수분해 문제와 연관시키는 것입니다. 알고리즘의 마지막 부분은 순전히 고전적인 부분으로, 양자 컴퓨터에서 위상 측정을 얻은 후 고전적인 컴퓨터에서 해결할 수 있습니다. 따라서 알고리즘의 마지막 부분은 실제 하드웨어에서 순서 찾기 회로를 실행하는 방법을 시연할 때까지 연기합니다.

하드웨어 실행

이제 이전에 ibm_marrakesh 에 대해 트랜스파일한 주문 찾기 회로를 실행할 수 있습니다. 여기서는 오류 억제를 위해 동적 디커플링 (DD)을, 오류 완화를 위해 게이트 트월링을 사용합니다. DD는 정확한 타이밍의 제어 펄스 시퀀스를 양자 디바이스에 적용하여 원치 않는 환경 상호작용과 비연속성을 효과적으로 평균화합니다. 반면 게이트 트월링은 특정 양자 게이트를 무작위화하여 코히어런트 오류를 폴리 오류로 변환하고, 이 오류는 2진법이 아닌 선형적으로 누적됩니다. 양자 계산의 일관성과 충실도를 높이기 위해 두 기술을 결합하는 경우가 많습니다.

# Sampler primitive to obtain the probability distribution
sampler = Sampler(backend)

# Turn on dynamical decoupling with sequence XpXm
sampler.options.dynamical_decoupling.enable = True
sampler.options.dynamical_decoupling.sequence_type = "XpXm"
# Enable gate twirling
sampler.options.twirling.enable_gates = True

pub = transpiled_circuit
job = sampler.run([pub], shots=1024)

result = job.result()[0]
counts = result.data["out"].get_counts()

plot_histogram(counts, figsize=(35, 5))

Output:

보시다시피, 가장 높은 개수의 동일한 비트스트링을 얻었습니다. 양자 하드웨어에는 노이즈가 있기 때문에 다른 비트스트링으로 일부 누수가 발생하는데, 이를 통계적으로 걸러낼 수 있습니다.

# Dictionary of bitstrings and their counts to keep
counts_keep = {}
# Threshold to filter
threshold = np.max(list(counts.values())) / 2

for key, value in counts.items():
    if value > threshold:
        counts_keep[key] = value

print(counts_keep)

Output:

{'00000000': 58, '01000000': 41, '11000000': 42, '10000000': 40}

4단계: 후처리 수행 및 원하는 클래식 형식으로 결과 반환

정수 인수분해

지금까지 위상 추정 회로를 사용하여 순서 찾기 문제를 구현하는 방법에 대해 설명했습니다. 이제 순서 찾기 문제를 정수 인수분해에 연결하여 쇼의 알고리즘을 완성합니다. 알고리즘의 이 부분은 고전적인 방식이라는 점에 유의하세요.

이제 $N = 15$ 와 $a = 2$ 의 예시를 통해 이를 증명해 보겠습니다. 우리가 측정한 위상은 $k / r$ 이며, $a^r \; (\textrm{mod} \; N) = 1$ 과 $k$ 는 $0$ 과 $r - 1$ 사이의 임의의 정수입니다. 이 방정식에서 $(a^r - 1) \; (\textrm{mod} \; N) = 0,$ 는 $N$ 이 $a^r-1$ 를 나누어야 함을 의미합니다. $r$ 도 짝수이면 $a^r -1 = (a^{r/2}-1)(a^{r/2}+1).$ 을 쓸 수 있습니다 $r$ 이 짝수가 아니면 더 나아갈 수 없고 $a$ 에 대해 다른 값으로 다시 시도해야 합니다; 그렇지 않으면 $N$ 과 $a^{r/2}-1$ 또는 $a^{r/2}+1$ 의 최대 공약수가 $N$ 의 고유 인수일 확률이 높습니다.

알고리즘의 일부 실행은 통계적으로 실패할 수 있으므로 $N$ 의 인자가 하나 이상 발견될 때까지 이 알고리즘을 반복합니다.

아래 셀은 $N=15$ 인자가 하나 이상 발견될 때까지 알고리즘을 반복합니다. 위의 하드웨어 실행 결과를 사용하여 각 반복에서 위상과 해당 요소를 추측합니다.

a = 2
N = 15

FACTOR_FOUND = False
num_attempt = 0

while not FACTOR_FOUND:
    print(f"\nATTEMPT {num_attempt}:")
    # Here, we get the bitstring by iterating over outcomes
    # of a previous hardware run with multiple shots.
    # Instead, we can also perform a single-shot measurement
    # here in the loop.
    bitstring = list(counts_keep.keys())[num_attempt]
    num_attempt += 1
    # Find the phase from measurement
    decimal = int(bitstring, 2)
    phase = decimal / (2**num_control)  # phase = k / r
    print(f"Phase: theta = {phase}")

    # Guess the order from phase
    frac = Fraction(phase).limit_denominator(N)
    r = frac.denominator  # order = r
    print(f"Order of {a} modulo {N} estimated as: r = {r}")

    if phase != 0:
        # Guesses for factors are gcd(a^{r / 2} ± 1, 15)
        if r % 2 == 0:
            x = pow(a, r // 2, N) - 1
            d = gcd(x, N)
            if d > 1:
                FACTOR_FOUND = True
                print(f"*** Non-trivial factor found: {x} ***")

Output:


ATTEMPT 0:
Phase: theta = 0.0
Order of 2 modulo 15 estimated as: r = 1

ATTEMPT 1:
Phase: theta = 0.25
Order of 2 modulo 15 estimated as: r = 4
*** Non-trivial factor found: 3 ***

토론

알고리즘의 확장성

RSA 암호화는 일반적으로 2048~4096비트 정도의 키 크기를 사용합니다. 쇼의 알고리즘으로 2048비트 숫자를 인수분해하려고 하면 오류 수정 오버헤드와 10억에 달하는 회로 깊이를 포함해 수백만 큐비트의 양자 회로가 생성되며, 이는 현재 양자 하드웨어의 실행 한계를 뛰어넘는 것입니다. 따라서 쇼의 알고리즘이 최신 암호화 시스템을 뚫기 위해서는 최적화된 회로 구성 방법이나 강력한 양자 오류 수정이 실질적으로 실행 가능해야 합니다. Shor 알고리즘의 리소스 추정에 대한 자세한 내용은 [9] 를 참조하세요.

문제점

튜토리얼을 완료하신 것을 축하드립니다! 지금이 바로 이해도를 테스트할 수 있는 좋은 시기입니다. 21을 인수분해하는 회로를 구성해 보시겠습니까? 원하는 $a$ 을 선택할 수 있습니다. 큐비트 수를 선택하기 위해 알고리즘의 비트 정확도를 결정해야 하며, 모듈식 지수 연산자 $M_a$ 를 설계해야 합니다. 이를 직접 시도해 본 다음 [6] 의 그림 9와 [7] 의 그림 2에 표시된 방법론을 읽어보시기 바랍니다.

def M_a_mod21():
    """
    M_a (mod 21)
    """

    # Your code here
    pass

참조

Shor, Peter W. "양자 컴퓨터에서 소인수분해 및 이산 로그에 대한 다항식 시간 알고리즘." SIAM 리뷰 41.2 (1999): 303-332.
IBM 존 와트루스 박사의 양자 알고리즘의 기초 강좌.
Vandersypen, Lieven MK, et al. " 핵자기공명을 이용한 쇼어의 양자소인수분해 알고리즘의 실험적 구현." 자연 414.6866 (2001): 883-887.
마틴-로페즈, 엔리케 외. " 큐비트 재활용을 이용한 쇼어의 양자 인수분해 알고리즘의 실험적 구현." 네이처 포토닉스 6.11 (2012): 773-776.
몬츠, 토마스 외. "확장 가능한 쇼 알고리즘의 실현." 과학 351.6277 (2016): 1068-1070.
아미코, 미르코, 자인 H. 살렘, 뮤어 쿰프. " IBM Q 경험을 사용한 Shor의 인수 분해 알고리즘에 대한 실험적 연구." 물리적 검토 A 100.1 (2019): 012305.
스코사나, 우나티, 마크 타메. " IBM 양자 프로세서에서 N=21 에 대한 Shor의 인수분해 알고리즘 데모." 과학 보고서 11.1 (2021): 16599.
키타에프, A. 유. "양자 측정과 아벨리안 안정기 문제." arXiv 사전 인쇄 quant-ph/9511026 (1995).
기드니, 크레이그, 마틴 에케라. "2천만 개의 노이즈 큐비트를 사용하여 8시간 안에 2048비트 RSA 정수를 인수분해하는 방법." 퀀텀 5 (2021): 433.

이 페이지가 도움이 되었습니까?

GitHub 에서 버그, 오타 또는 콘텐츠 요청을 신고하세요.

쇼어의 알고리즘

요구사항

설정

1단계: 고전적 입력을 양자 문제에 매핑하기

배경

위상 추정 문제

주문 찾기 문제

모듈러 지수 연산자

2단계: 양자 하드웨어 실행을 위한 문제 최적화

N=15N = 15N=15 및 a=2a=2a=2 의 구체적인 예시

3단계: Qiskit primitives 명령어로 실행합니다

하드웨어 실행

4단계: 후처리 수행 및 원하는 클래식 형식으로 결과 반환

정수 인수분해

토론

관련 작업

알고리즘의 확장성

문제점

참조

$N = 15$ 및 $a=2$ 의 구체적인 예시

3단계: `Qiskit primitives` 명령어로 실행합니다