Algoritmo ni Shor

Tansiya ng paggamit: Tatlong segundo sa isang Eagle r3 processor (TANDAAN: Ito ay isang tantiya lamang. Maaaring mag-iba ang inyong runtime.)

Ang algoritmo ni Shor, na binuo ni Peter Shor noong 1994, ay isang makabagong quantum algorithm para sa pag-factor ng mga integer sa polynomial time. Ang kahalagahan nito ay nakasalalay sa kakayahan nitong mag-factor ng malalaking integer nang mas mabilis nang exponential kaysa sa alinmang kilalang classical algorithm, na bumabanta sa seguridad ng malawakang ginagamit na mga cryptographic system tulad ng RSA, na umaasa sa kahirapan ng pag-factor ng malalaking numero. Sa pamamagitan ng epektibong paglutas ng problemang ito sa sapat na makapangyarihang quantum computer, ang algoritmo ni Shor ay maaaring magdulot ng rebolusyon sa mga larangan tulad ng cryptography, cybersecurity, at computational mathematics, na nagpapahayag ng transformative power ng quantum computation.

Nakatuon ang tutorial na ito sa pagpapakita ng algoritmo ni Shor sa pamamagitan ng pag-factor ng 15 sa isang quantum computer.

Una, tutukuyin natin ang order finding problem at bubuo ng mga katumbas na circuit mula sa quantum phase estimation protocol. Susunod, pagaganahin natin ang mga order finding circuit sa tunay na hardware gamit ang pinakamaikling-lalim na mga circuit na maaari nating i-transpile. Ang huling seksyon ay kukumpletuhin ang algoritmo ni Shor sa pamamagitan ng pag-uugnay ng order finding problem sa integer factorization.

Tatapusin natin ang tutorial na may pag-uusap tungkol sa iba pang mga demonstrasyon ng algoritmo ni Shor sa tunay na hardware, na nakatuon sa parehong generic na mga implementasyon at sa mga iniayon para sa pag-factor ng mga partikular na integer tulad ng 15 at 21. Tandaan: Ang tutorial na ito ay mas nakatuon sa implementasyon at demonstrasyon ng mga circuit na may kaugnayan sa algoritmo ni Shor. Para sa isang malalim na educational resource sa materyal, mangyaring sumangguni sa kursong Fundamentals of quantum algorithms ni Dr. John Watrous, at mga papel sa seksyong References.

Mga Kinakailangan

Bago magsimula sa tutorial na ito, siguraduhin na mayroon kayo ng mga sumusunod na naka-install:

• Qiskit SDK v2.0 o mas bago, na may suportang visualization
• Qiskit Runtime v0.40 o mas bago (pip install qiskit-ibm-runtime)

Setup

# Added by doQumentation — required packages for this notebook
!pip install -q numpy pandas qiskit 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

Hakbang 1: I-map ang mga classical input sa isang quantum problem

Background

Ang algoritmo ni Shor para sa integer factorization ay gumagamit ng isang intermediary problem na kilala bilang order finding problem. Sa seksyong ito, ipapakita natin kung paano lutasin ang order finding problem gamit ang quantum phase estimation.

Problema sa phase estimation

Sa phase estimation problem, binibigyan tayo ng isang quantum state na $\ket{\psi}$ ng $n$ qubit, kasama ng isang unitary quantum circuit na gumagana sa $n$ qubit. Ipinangako sa atin na ang $\ket{\psi}$ ay isang eigenvector ng unitary matrix na $U$ na naglalarawan sa aksyon ng circuit, at ang ating layunin ay kalkulahin o tantiyahin ang eigenvalue na $\lambda = e^{2 \pi i \theta}$ kung saan tumutugma ang $\ket{\psi}$. Sa madaling salita, ang circuit ay dapat mag-output ng approximation sa numerong $\theta \in [0, 1)$ na nakakatugon sa $U \ket{\psi}= e^{2 \pi i \theta} \ket{\psi}.$ Ang layunin ng phase estimation circuit ay tantiyahin ang $\theta$ sa $m$ bits. Sa matematika, nais nating makahanap ng $y$ na $\theta \approx y / 2^m$, kung saan ang $y \in {0, 1, 2, \dots, 2^{m-1}}$. Ang sumusunod na larawan ay nagpapakita ng quantum circuit na tumantiya ng $y$ sa $m$ bits sa pamamagitan ng paggawa ng measurement sa $m$ qubit. Sa circuit sa itaas, ang itaas na $m$ qubit ay sinimulan sa estado ng $\ket{0^m}$, at ang ibabang $n$ qubit ay sinimulan sa $\ket{\psi}$, na ipinangakong isang eigenvector ng $U$. Ang unang sangkap sa phase estimation circuit ay ang mga controlled-unitary operation na responsable sa pagsasagawa ng phase kickback sa kanilang katumbas na control qubit. Ang mga controlled unitary na ito ay naka-exponentiate alinsunod sa posisyon ng control qubit, mula sa least significant bit hanggang sa most significant bit. Dahil ang $\ket{\psi}$ ay isang eigenvector ng $U$, ang estado ng ibabang $n$ qubit ay hindi naapektuhan ng operasyong ito, ngunit ang phase information ng eigenvalue ay kumakalat sa itaas na $m$ qubit. Lumalabas na pagkatapos ng phase kickback operation sa pamamagitan ng mga controlled-unitary, ang lahat ng posibleng estado ng itaas na $m$ qubit ay orthonormal sa isa't isa para sa bawat eigenvector na $\ket{\psi}$ ng unitary na $U$. Samakatuwid, ang mga estadong ito ay perpektong maipag-kakaiba, at maaari nating i-rotate ang basis na kanilang binubuo pabalik sa computational basis upang gumawa ng measurement. Ang isang mathematical analysis ay nagpapakita na ang rotation matrix na ito ay tumutugma sa inverse quantum Fourier transform (QFT) sa $2^m$-dimensional Hilbert space. Ang intuition sa likod nito ay ang periodic structure ng mga modular exponentiation operator ay naka-encode sa quantum state, at ang QFT ay nag-convert ng periodicity na ito sa mga nasu-sukat na peak sa frequency domain.

Para sa mas malalim na pag-unawa kung bakit ginagamit ang QFT circuit sa algoritmo ni Shor, inirerekomenda namin sa mambabasa ang kursong Fundamentals of quantum algorithms. Handa na tayong gumamit ng phase estimation circuit para sa order finding.

Problema sa order finding

Upang tukuyin ang order finding problem, magsisimula tayo sa ilang konsepto ng number theory. Una, para sa anumang ibinigay na positive integer na $N$, tukuyin ang set na $\mathbb{Z}_N$ bilang $\mathbb{Z}_N = \{0, 1, 2, \dots, N-1\}.$ Ang lahat ng arithmetic operation sa $\mathbb{Z}_N$ ay ginagawa modulo $N$. Partikular, ang lahat ng elemento na $a \in \mathbb{Z}_n$ na coprime sa $N$ ay espesyal at bumubuo ng $\mathbb{Z}^*_N$ bilang $\mathbb{Z}^*_N = \{ a \in \mathbb{Z}_N : \mathrm{gcd}(a, N)=1 \}.$ Para sa isang elemento na $a \in \mathbb{Z}^*_N$, ang pinakamaliit na positive integer na $r$ na $a^r \equiv 1 \; (\mathrm{mod} \; N)$ ay tinutukoy bilang order ng $a$ modulo $N$. Tulad ng makikita natin sa ibang pagkakataon, ang paghahanap ng order ng isang $a \in \mathbb{Z}^*_N$ ay magbibigay-daan sa atin upang i-factor ang $N$. Upang buuin ang order finding circuit mula sa phase estimation circuit, kailangan natin ng dalawang pagsasaalang-alang. Una, kailangan nating tukuyin ang unitary na $U$ na magpapahintulot sa atin na mahanap ang order na $r$, at pangalawa, kailangan nating tukuyin ang isang eigenvector na $\ket{\psi}$ ng $U$ upang ihanda ang panimulang estado ng phase estimation circuit.

Upang iugnay ang order finding problem sa phase estimation, isasaalang-alang natin ang operasyon na tinukoy sa isang system na ang mga classical state ay tumutugma sa $\mathbb{Z}_N$, kung saan tayo ay nag-multiply ng isang fixed element na $a \in \mathbb{Z}^*_N$. Partikular, tutukuyin natin ang multiplication operator na $M_a$ na $M_a \ket{x} = \ket{ax \; (\mathrm{mod} \; N)}$ para sa bawat $x \in \mathbb{Z}_N$. Tandaan na implisito na ang produkto modulo $N$ ay kinukuha natin sa loob ng ket sa kanang bahagi ng equation. Ang isang mathematical analysis ay nagpapakita na ang $M_a$ ay isang unitary operator. Higit pa rito, lumalabas na ang $M_a$ ay may mga pares ng eigenvector at eigenvalue na nagpapahintulot sa atin na iugnay ang order na $r$ ng $a$ sa phase estimation problem. Partikular, para sa anumang pagpili ng $j \in \{0, \dots, r-1\}$, mayroon tayo na $\ket{\psi_j} = \frac{1}{\sqrt{r}} \sum^{r-1}_{k=0} \omega^{-jk}_{r} \ket{a^k}$ ay isang eigenvector ng $M_a$ na ang katumbas na eigenvalue ay $\omega^{j}_{r}$, kung saan $\omega^{j}_{r} = e^{2 \pi i \frac{j}{r}}.$ Sa pag-obserba, nakikita natin na ang isang convenient na pares ng eigenvector/eigenvalue ay ang estado na $\ket{\psi_1}$ na may $\omega^{1}_{r} = e^{2 \pi i \frac{1}{r}}$. Samakatuwid, kung makakahanap tayo ng eigenvector na $\ket{\psi_1}$, maaari nating tantiyahin ang phase na $\theta=1/r$ gamit ang ating quantum circuit at sa gayon ay makakuha ng tantiya ng order na $r$. Gayunpaman, hindi ito madali gawin, at kailangan nating isaalang-alang ang alternatibo.

Isaalang-alang natin kung ano ang magreresulta sa circuit kung ihanda natin ang computational state na $\ket{1}$ bilang panimulang estado. Ito ay hindi isang eigenstate ng $M_a$, ngunit ito ay ang uniform superposition ng mga eigenstate na aming inilarawan sa itaas. Sa madaling salita, ang sumusunod na relasyon ay totoo. $\ket{1} = \frac{1}{\sqrt{r}} \sum^{r-1}_{k=0} \ket{\psi_k}$ Ang implikasyon ng equation sa itaas ay kung itakda natin ang panimulang estado sa $\ket{1}$, makakakuha tayo ng eksaktong parehong resulta ng measurement na parang pumili tayo ng $k \in \{ 0, \dots, r-1\}$ nang pantay-pantay na random at ginamit ang $\ket{\psi_k}$ bilang eigenvector sa phase estimation circuit. Sa madaling salita, ang isang measurement ng itaas na $m$ qubit ay nagbubunga ng approximation na $y / 2^m$ sa halaga ng $k / r$ kung saan ang $k \in \{ 0, \dots, r-1\}$ ay pinili nang pantay-pantay na random. Ito ay nagbibigay-daan sa atin na matutunan ang $r$ na may mataas na antas ng kumpiyansa pagkatapos ng ilang independyenteng run, na siyang ating layunin.

Mga modular exponentiation operator

Sa ngayon, inugnay natin ang phase estimation problem sa order finding problem sa pamamagitan ng pagtukoy ng $U = M_a$ at $\ket{\psi} = \ket{1}$ sa ating quantum circuit. Samakatuwid, ang huling natitirang sangkap ay ang paghahanap ng isang epektibong paraan upang tukuyin ang mga modular exponential ng $M_a$ bilang $M_a^k$ para sa $k = 1, 2, 4, \dots, 2^{m-1}$. Upang isagawa ang computation na ito, nakikita natin na para sa anumang power na $k$ na ating pipiliin, maaari tayong lumikha ng circuit para sa $M_a^k$ hindi sa pamamagitan ng pag-uulit ng $k$ beses ang circuit para sa $M_a$, kundi sa halip sa pamamagitan ng pagkalkula ng $b = a^k \; \mathrm{mod} \; N$ at pagkatapos ay paggamit ng circuit para sa $M_b$. Dahil kailangan lang natin ng mga power na mga power din ng 2, maaari nating gawin ito nang epektibo sa classical gamit ang iterative squaring.

Hakbang 2: I-optimize ang problema para sa quantum hardware execution

Partikular na halimbawa na may $N = 15$ at $a=2$

Maaari tayong huminto dito upang talakayin ang isang partikular na halimbawa at buuin ang order finding circuit para sa $N=15$. Tandaan na ang posibleng nontrivial na $a \in \mathbb{Z}_N^*$ para sa $N=15$ ay $a \in \{2, 4, 7, 8, 11, 13, 14 \}$. Para sa halimbawang ito, pipiliin natin ang $a=2$. Bubuo tayo ng $M_2$ operator at ng mga modular exponentiation operator na $M_2^k$. Ang aksyon ng $M_2$ sa mga computational basis state ay ang mga sumusunod. $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}$ Sa pag-obserba, nakikita natin na ang mga basis state ay na-shuffle, kaya mayroon tayong permutation matrix. Maaari nating buuin ang operasyong ito sa apat na qubit gamit ang mga swap gate. Sa ibaba, bubuo tayo ng $M_2$ at ng mga controlled-$M_2$ operation.

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)

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)

Ang mga gate na gumagana sa mahigit dalawang qubit ay pag-hahatiin pa sa mga two-qubit gate.

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

Ngayon kailangan nating buuin ang mga modular exponentiation operator. Upang makakuha ng sapat na precision sa phase estimation, gagamitin natin ang walong qubit para sa estimation measurement. Samakatuwid, kailangan nating buuin ang $M_b$ na may $b = a^{2^k} \; (\mathrm{mod} \; N)$ para sa bawat $k = 0, 1, \dots, 7$.

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)

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

Tulad ng makikita natin mula sa listahan ng mga halaga ng $b$, bilang karagdagan sa $M_2$ na dati nating ginawa, kailangan din nating bumuo ng $M_4$ at $M_1$. Tandaan na ang $M_1$ ay gumagana nang trivial sa mga computational basis state, kaya ito ay simpleng identity operator.

Ang $M_4$ ay gumagana sa mga computational basis state na ganito. $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}$

Samakatuwid, ang permutation na ito ay maaaring buuin gamit ang sumusunod na swap operation.

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)

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)

Ang mga gate na gumagana sa mahigit dalawang qubit ay pag-hahatiin pa sa mga two-qubit gate.

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

Nakita natin na ang mga $M_b$ operator para sa isang ibinigay na $b \in \mathbb{Z}^*_N$ ay mga permutation operation. Dahil sa relatibong maliit na laki ng permutation problem na mayroon tayo dito, dahil ang $N=15$ ay nangangailangan lamang ng apat na qubit, nakagawa tayo ng mga operasyong ito nang direkta gamit ang mga SWAP gate sa pamamagitan ng inspeksyon. Sa pangkalahatan, maaaring hindi ito scalable na pamamaraan. Sa halip, maaaring kailangan nating buuin nang hayagan ang permutation matrix, at gamitin ang UnitaryGate class ng Qiskit at mga transpilation method upang i-synthesize ang permutation matrix na ito. Gayunpaman, maaari itong magresulta sa mas malalim na circuit. Ang isang halimbawa ay sumusunod.

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
)

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})

Ihambing natin ang mga bilang na ito sa compiled circuit depth ng ating manual na implementasyon ng $M_2$ gate.

# 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
)

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

Tulad ng makikita natin, ang permutation matrix approach ay nagresulta sa mas malalim na circuit kahit para sa isang $M_2$ gate lamang kumpara sa ating manual na implementasyon nito. Samakatuwid, magpapatuloy tayo sa ating naunang implementasyon ng mga $M_b$ operation. Ngayon, handa na tayong buuin ang buong order finding circuit gamit ang ating mga dati nang tinukoy na controlled modular exponentiation operator. Sa sumusunod na code, ini-import din natin ang QFT circuit mula sa Qiskit Circuit library, na gumagamit ng mga Hadamard gate sa bawat qubit, isang serye ng mga controlled-U1 (o Z, depende sa phase) gate, at isang layer ng mga swap gate.

# 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)

Tandaan na inalis natin ang mga controlled modular exponentiation operation mula sa natitirang control qubit dahil ang $M_1$ ay identity operator. Tandaan na sa ibang pagkakataon sa tutorial na ito, pagaganahin natin ang circuit na ito sa ibm_marrakish backend. Upang gawin ito, i-transpile natin ang circuit ayon sa partikular na backend na ito at iulat ang circuit depth at mga gate count.

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
)

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

Hakbang 3: Magsagawa gamit ang mga Qiskit primitive

Una, talakayin natin kung ano ang teoretikal na makukuha natin kung pagaganahin natin ang circuit na ito sa isang ideal na simulator. Sa ibaba, mayroon tayong set ng mga resulta ng simulation ng circuit sa itaas gamit ang 1024 shot. Tulad ng makikita natin, nakakakuha tayo ng halos pantay na pamamahagi sa apat na bitstring sa mga control qubit.

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

plot_histogram(counts)

Sa pamamagitan ng pagsukatan sa mga control qubit, nakakakuha tayo ng walong-bit na phase estimation ng $M_a$ operator. Maaari nating i-convert ang binary representation na ito sa decimal upang mahanap ang nasukat na phase. Tulad ng makikita natin mula sa histogram sa itaas, apat na iba't ibang bitstring ang nasukat, at ang bawat isa sa kanila ay tumutugma sa isang phase value na ganito.

# 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)

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

Alalahanin na ang anumang nasukat na phase ay tumutugma sa $\theta = k / r$ kung saan ang $k$ ay sample nang pantay-pantay na random mula sa $\{0, 1, \dots, r-1 \}$. Samakatuwid, maaari nating gamitin ang continued fractions algorithm upang subukang mahanap ang $k$ at ang order na $r$. Ang Python ay may built-in na functionality na ito. Maaari nating gamitin ang fractions module upang gawing Fraction object ang isang float, halimbawa:

Fraction(0.666)

Fraction(5998794703657501, 9007199254740992)

Dahil nagbibigay ito ng mga fraction na ibinabalik ang resulta nang eksakto (sa kasong ito, 0.6660000...), maaari itong magbigay ng mga komplikadong resulta tulad ng nasa itaas. Maaari nating gamitin ang .limit_denominator() method upang makuha ang fraction na pinaka-katulad ng ating float, na may denominator sa ibaba ng isang tiyak na halaga:

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

Fraction(2, 3)

Ito ay mas maganda. Ang order (r) ay dapat na mas mababa sa N, kaya itakda natin ang maximum denominator na 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)

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

Makikita natin na dalawa sa mga nasukat na eigenvalue ay nagbigay sa atin ng tamang resulta: $r=4$, at makikita natin na ang algoritmo ni Shor para sa order finding ay may tsansang mabigo. Ang mga masamang resulta na ito ay dahil sa $k = 0$, o dahil ang $k$ at $r$ ay hindi coprime - at sa halip na $r$, binibigyan tayo ng isang factor ng $r$. Ang pinakamadaling solusyon dito ay simpleng ulitin ang eksperimento hanggang sa makakuha tayo ng kasiya-siyang resulta para sa $r$. Sa ngayon, naipatupad na natin ang order finding problem para sa $N=15$ na may $a=2$ gamit ang phase estimation circuit sa isang simulator. Ang huling hakbang ng algoritmo ni Shor ay ang pag-ugnay ng order finding problem sa integer factorization problem. Ang huling bahaging ito ng algorithm ay purong classical at maaaring lutasin sa isang classical computer pagkatapos makuha ang mga phase measurement mula sa isang quantum computer. Samakatuwid, ipinagpaliban natin ang huling bahagi ng algorithm hanggang sa ipakita natin kung paano natin maaaring paganahin ang order finding circuit sa tunay na hardware.

Mga hardware run

Ngayon ay maaari nating paganahin ang order finding circuit na dati nating na-transpile para sa ibm_marrakish. Dito ginagamit natin ang dynamical decoupling (DD) para sa error suppression, at gate twirling para sa mga layunin ng error mitigation. Ang DD ay nagsasangkot ng pag-apply ng mga sequence ng tumpak na naka-time na control pulse sa isang quantum device, na epektibong nag-aaverage out ng mga hindi gustong environmental interaction at decoherence. Ang gate twirling, sa kabilang banda, ay nag-randomize ng mga partikular na quantum gate upang gawing Pauli error ang mga coherent error, na nag-iipon nang linear sa halip na quadratic. Ang parehong pamamaraan ay madalas na pinagsama upang pahusayin ang coherence at fidelity ng mga quantum computation.

# 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))

Tulad ng makikita natin, nakuha natin ang parehong mga bitstring na may pinakamataas na bilang. Dahil ang quantum hardware ay may ingay, may ilang leakage sa ibang mga bitstring, na maaari nating i-filter sa pamamagitan ng estadistika.

# 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)

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

Hakbang 4: Mag-post-process at ibalik ang resulta sa ninanais na classical format

Pag-factor ng Integer

Sa ngayon, tinalakay natin kung paano natin maaaring ipatupad ang order finding problem gamit ang isang phase estimation circuit. Ngayon, iuugnay natin ang order finding problem sa integer factorization, na kumukumpleto sa algoritmo ni Shor. Tandaan na ang bahaging ito ng algorithm ay classical. Ipakikita natin ito gamit ang ating halimbawa ng $N = 15$ at $a = 2$. Alalahanin na ang phase na ating nasukat ay $k / r$, kung saan $a^r \; (\textrm{mod} \; N) = 1$ at ang $k$ ay isang random integer sa pagitan ng $0$ at $r - 1$. Mula sa equation na ito, mayroon tayo $(a^r - 1) \; (\textrm{mod} \; N) = 0,$ na nangangahulugan na ang $N$ ay dapat na naghahati sa $a^r-1$. Kung ang $r$ ay even din, kung gayon maaari nating isulat $a^r -1 = (a^{r/2}-1)(a^{r/2}+1).$ Kung ang $r$ ay hindi even, hindi na tayo makakapag-patuloy at kailangang subukan muli na may ibang halaga para sa $a$; kung hindi man, may mataas na probabilidad na ang greatest common divisor ng $N$ at alinman sa $a^{r/2}-1$, o $a^{r/2}+1$ ay isang wastong factor ng $N$.

Dahil ang ilang run ng algorithm ay mabibigo sa estadistika, uulitin natin ang algorithm na ito hanggang sa makahanap ng hindi bababa sa isang factor ng $N$. Ang cell sa ibaba ay umuulit ng algorithm hanggang sa makahanap ng hindi bababa sa isang factor ng $N=15$. Gagamitin natin ang mga resulta ng hardware run sa itaas upang tantiyahin ang phase at ang katumbas na factor sa bawat iteration.

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} ***")

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 ***

Talakayan

Kaugnay na gawain

Sa seksyong ito, tinalakay natin ang iba pang milestone na gawain na nagpakita ng algoritmo ni Shor sa tunay na hardware.

Ang seminal work na [3] mula sa IBM® ay nagpakita ng algoritmo ni Shor sa unang pagkakataon, na nag-factor ng numerong 15 sa mga prime factor nitong 3 at 5 gamit ang isang seven-qubit nuclear magnetic resonance (NMR) quantum computer. Ang isa pang eksperimento [4] ay nag-factor ng 15 gamit ang mga photonic qubit. Sa pamamagitan ng paggamit ng isang qubit na na-recycle nang maraming beses at pag-encode ng work register sa mas mataas-dimensional na estado, binawasan ng mga mananaliksik ang kinakailangang bilang ng qubit sa isang-katlo ng nasa standard protocol, gamit ang isang two-photon compiled algorithm. Ang isang makabuluhang papel sa demonstrasyon ng algoritmo ni Shor ay [5], na gumagamit ng iterative phase estimation technique ni Kitaev [8] upang mabawasan ang qubit requirement ng algorithm. Gumamit ang mga may-akda ng pitong control qubit at apat na cache qubit, kasama ang implementasyon ng mga modular multiplier. Ang implementasyong ito, gayunpaman, ay nangangailangan ng mid-circuit measurement na may feed-forward operation at qubit recycling na may reset operation. Ang demonstrasyong ito ay ginawa sa isang ion-trap quantum computer.

Ang mas kamakailang gawain [6] ay nakatuon sa pag-factor ng 15, 21, at 35 sa IBM Quantum® hardware. Katulad ng nakaraang gawain, gumamit ang mga mananaliksik ng compiled version ng algorithm na gumamit ng semi-classical quantum Fourier transform gaya ng iminungkahi ni Kitaev upang mabawasan ang bilang ng mga physical qubit at gate. Ang pinakabagong gawain [7] ay nagsagawa din ng proof-of-concept demonstration para sa pag-factor ng integer na 21. Ang demonstrasyong ito ay nagsangkot din ng paggamit ng compiled version ng quantum phase estimation routine, at itinayo ang nakaraang demonstrasyon sa pamamagitan ng [4]. Lumampas ang mga may-akda sa gawain na ito sa pamamagitan ng paggamit ng configuration ng mga approximate Toffoli gate na may residual phase shift. Ang algorithm ay na-implement sa mga IBM quantum processor gamit lamang ang limang qubit, at ang presensya ng entanglement sa pagitan ng control at register qubit ay matagumpay na napatunayan.

Scaling ng algorithm

Tandaan natin na ang RSA encryption ay karaniwang nagsasangkot ng mga key size na nasa order ng 2048 hanggang 4096 bit. Ang pagtatangka na i-factor ang isang 2048-bit na numero gamit ang algoritmo ni Shor ay magreresulta sa isang quantum circuit na may milyun-milyong qubit, kabilang ang error correction overhead at isang circuit depth na nasa order ng isang bilyon, na lampas sa mga limitasyon ng kasalukuyang quantum hardware na maisagawa. Samakatuwid, ang algoritmo ni Shor ay mangangailangan ng alinman sa mga na-optimize na circuit construction method o matatag na quantum error correction upang maging praktikal na maabot para sa pagwasak ng mga modernong cryptographic system. Inirerekomenda namin kayo sa [9] para sa mas detalyadong talakayan sa resource estimation para sa algoritmo ni Shor.

Hamon

Binabati namin kayo sa pagtapos ng tutorial! Ngayon ay mahusay na panahon upang subukan ang inyong pag-unawa. Maaari ba kayong subukang buuin ang circuit para sa pag-factor ng 21? Maaari kayong pumili ng sariling $a$. Kakailanganin ninyong magpasya sa bit accuracy ng algorithm upang piliin ang bilang ng mga qubit, at kakailanganin ninyong idisenyo ang mga modular exponentiation operator na $M_a$. Hinihikayat namin kayong subukan ito mismo, at pagkatapos ay basahin ang tungkol sa mga metodolohiya na ipinakita sa Fig. 9 ng [6] at Fig. 2 ng [7].

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

    # Your code here
    pass

Mga Reperensya

1. Shor, Peter W. "Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer." SIAM review 41.2 (1999): 303-332.
2. IBM Quantum Fundamentals of Quantum Algorithms course by Dr. John Watrous.
3. Vandersypen, Lieven MK, et al. "Experimental realization of Shor's quantum factoring algorithm using nuclear magnetic resonance." Nature 414.6866 (2001): 883-887.
4. Martin-Lopez, Enrique, et al. "Experimental realization of Shor's quantum factoring algorithm using qubit recycling." Nature photonics 6.11 (2012): 773-776.
5. Monz, Thomas, et al. "Realization of a scalable Shor algorithm." Science 351.6277 (2016): 1068-1070.
6. Amico, Mirko, Zain H. Saleem, and Muir Kumph. "Experimental study of Shor's factoring algorithm using the IBM Q Experience." Physical Review A 100.1 (2019): 012305.
7. Skosana, Unathi, and Mark Tame. "Demonstration of Shor's factoring algorithm for N=21 on IBM quantum processors." Scientific reports 11.1 (2021): 16599.
8. Kitaev, A. Yu. "Quantum measurements and the Abelian stabilizer problem." arXiv preprint quant-ph/9511026 (1995).
9. Gidney, Craig, and Martin Ekerå. "How to factor 2048 bit RSA integers in 8 hours using 20 million noisy qubits." Quantum 5 (2021): 433.