Algorytmy kwantowe: Estymacja fazy

uwaga

Kento Ueda (15 maja 2024)

Ten notatnik przedstawia podstawowe koncepcje oraz implementację Kwantowej Transformacji Fouriera (QFT) i Kwantowej Estymacji Fazy (QPE).

Pobierz pdf oryginalnego wykładu. Uwaga: niektóre fragmenty kodu mogą być przestarzałe, ponieważ są to statyczne obrazy.

Przybliżony czas QPU potrzebny do uruchomienia tego eksperymentu wynosi 7 sekund.

1. Wprowadzenie

Kwantowa Transformacja Fouriera (QFT)

Kwantowa Transformacja Fouriera jest kwantowym odpowiednikiem klasycznej dyskretnej transformacji Fouriera. Jest to transformacja liniowa stosowana do stanów kwantowych, odwzorowująca bazy obliczeniowe na ich reprezentacje w bazie Fouriera. QFT odgrywa kluczową rolę w wielu algorytmach kwantowych, oferując efektywną metodę wydobywania informacji o okresowości ze stanów kwantowych. QFT można zaimplementować z $O(N^2)$ operacjami przy użyciu bramek kwantowych takich jak bramki Hadamarda i bramki Control-Phase dla $N$ kubitów, co umożliwia wykładnicze przyspieszenie w porównaniu z klasyczną transformacją Fouriera.

• Zastosowania: Jest fundamentalnym elementem algorytmów kwantowych, takich jak algorytm Shora do faktoryzacji dużych liczb całkowitych oraz logarytmu dyskretnego.

Kwantowa Estymacja Fazy (QPE)

Kwantowa Estymacja Fazy to algorytm kwantowy wykorzystywany do estymacji fazy związanej z wektorem własnym operatora unitarnego. Algorytm ten stanowi pomost między abstrakcyjnymi właściwościami matematycznymi stanów kwantowych a ich zastosowaniami obliczeniowymi.

• Zastosowania: Może rozwiązywać problemy takie jak znajdowanie wartości własnych macierzy unitarnych i symulacja układów kwantowych.

Razem QFT i QPE stanowią niezbędny trzon wielu algorytmów kwantowych rozwiązujących problemy, które są niewykonalne dla komputerów klasycznych. Po przeczytaniu tego notatnika zrozumiesz, w jaki sposób te techniki są implementowane.

2. Podstawy Kwantowej Transformacji Fouriera (QFT)

# Added by doQumentation — required packages for this notebook
!pip install -q numpy qiskit qiskit-aer qiskit-ibm-runtime

from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister
from qiskit_aer import AerSimulator
from qiskit.visualization import plot_histogram, plot_bloch_multivector
from qiskit.quantum_info import Statevector
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit_ibm_runtime import Sampler

from numpy import pi

Przez analogię do dyskretnej transformacji Fouriera, QFT działa na stanie kwantowym $\vert X\rangle = \sum_{j=0}^{N-1} x_j \vert j \rangle$ dla $N$ kubitów i odwzorowuje go na stan kwantowy $\vert Y\rangle = \sum_{k=0}^{N-1} y_k \vert k \rangle$.

$$QFT_{N}: \vert j \rangle \mapsto \frac{1}{\sqrt{N}}\sum_{k=0}^{N-1}\omega_N^{jk} \vert k \rangle$$

gdzie $\omega_N^{jk} = e^{2\pi i \frac{jk}{N}}$.

Lub reprezentacja w postaci macierzy unitarnej:

$$U_{QFT} = \frac{1}{\sqrt{N}} \sum_{j=0}^{N-1} \sum_{k=0}^{N-1} \omega_N^{jk} \vert k \rangle \langle j \vert$$

2.1. Intuicja

Kwantowa transformacja Fouriera (QFT) dokonuje przekształcenia między dwiema bazami: bazą obliczeniową (Z) i bazą Fouriera. Ale co oznacza baza Fouriera w tym kontekście? Prawdopodobnie pamiętasz, że transformata Fouriera $F(\omega)$ funkcji $f(x)$ opisuje splot $f(x)$ z funkcją sinusoidalną o częstotliwości $\omega$. Mówiąc prosto: transformata Fouriera to funkcja opisująca, ile z każdej częstotliwości $\omega$ potrzebowalibyśmy, aby zbudować funkcję $f(x)$ z funkcji sinus (lub cosinus). Aby lepiej zrozumieć, co oznacza QFT w tym kontekście, rozważ poniższe obrazy krokowe pokazujące liczbę zakodowaną binarnie, używając czterech kubitów:

W bazie obliczeniowej przechowujemy liczby binarnie, używając stanów $|0\rangle$ i $|1\rangle$.

Zwróć uwagę na częstotliwość, z jaką zmieniają się poszczególne kubity; skrajnie lewy kubit zmienia się przy każdym zwiększeniu liczby, kolejny co 2 zwiększenia, trzeci co 4 zwiększenia i tak dalej.

Jeśli zastosujemy kwantową transformację Fouriera do tych stanów, odwzorowujemy

$$|\text{Stan w bazie obliczeniowej}\rangle \quad \xrightarrow[]{\text{QFT}} \quad |\text{Stan w bazie Fouriera}\rangle$$ $$\text{QFT}|x\rangle = |\widetilde{x}\rangle$$

(Często oznaczamy stany w bazie Fouriera za pomocą tyldy (~)).

W bazie Fouriera przechowujemy liczby, używając różnych obrotów wokół osi Z:

IFRAME

Liczba, którą chcemy przechować, określa kąt, pod jakim każdy kubit jest obracany wokół osi Z. W stanie $|\widetilde{0}\rangle$ wszystkie kubity znajdują się w stanie $|{+}\rangle$. Jak widać w powyższym przykładzie, aby zakodować stan $|\widetilde{5}\rangle$ na 4 kubitach, obróciliśmy skrajnie lewy kubit o $\tfrac{5}{2^n} = \tfrac{5}{16}$ pełnego obrotu ($\tfrac{5}{16}\times 2\pi$ radianów). Następny kubit jest obrócony dwukrotnie więcej ($\tfrac{10}{16}\times 2\pi$ radianów, czyli $10/16$ pełnego obrotu), następnie kąt ten jest podwajany dla kolejnego kubitu i tak dalej.

Ponownie, zwróć uwagę na częstotliwość, z jaką zmienia się każdy kubit. Skrajnie lewy kubit (qubit 0) w tym przypadku ma najniższą częstotliwość, a skrajnie prawy — najwyższą.

2.2 Przykład: 1-kubitowa QFT

Rozważmy przypadek $N = 2^1 = 2$.

Macierz unitarną można zapisać:

$$\begin{aligned} U_{QFT} & = \frac{1}{\sqrt{2}} \sum_{j=0}^{1} \sum_{k=0}^{1} \omega_2^{jk} \vert k \rangle \langle j \vert \\ & = \frac{1}{\sqrt{2}} (\omega_2^{0} \vert 0 \rangle \langle 0 \vert + \omega_2^{0} \vert 0 \rangle \langle 1 \vert + \omega_2^{0} \vert 1 \rangle \langle 0 \vert + \omega_2^{1} \vert 1 \rangle \langle 1 \vert) \\ & = \frac{1}{\sqrt{2}} (\vert 0 \rangle \langle 0 \vert + \vert 0 \rangle \langle 1 \vert + \vert 1 \rangle \langle 0 \vert - \vert 1 \rangle \langle 1 \vert) \\ & = H \end{aligned}$$

Ta operacja jest wynikiem zastosowania bramki Hadamarda ($H$).

2.3 Reprezentacja produktowa QFT

Uogólnijmy transformację dla $N = 2^n$, $QFT_{N}$ działającą na stanie $\vert x \rangle = \vert x_1\ldots x_n \rangle$.

$$\begin{aligned} QFT_N\vert x \rangle & = \frac{1}{\sqrt{N}} \sum_{y=0}^{N-1}\omega_N^{xy} \vert y \rangle \\ & = \frac{1}{\sqrt{N}} \sum_{y=0}^{N-1} e^{2 \pi i xy / 2^n} \vert y \rangle ~\text{ponieważ}\: \omega_N^{xy} = e^{2\pi i \frac{xy}{N}} \:\text{oraz}\: N = 2^n \\ & = \frac{1}{\sqrt{N}} \sum_{y=0}^{N-1} e^{2 \pi i \left(\sum_{k=1}^n y_k/2^k\right) x} \vert y_1 \ldots y_n \rangle \:\text{zapisując w ułamkowej notacji binarnej}\: y = y_1\ldots y_n, y/2^n = \sum_{k=1}^n y_k/2^k \\ & = \frac{1}{\sqrt{N}} \sum_{y=0}^{N-1} \prod_{k=1}^n e^{2 \pi i x y_k/2^k } \vert y_1 \ldots y_n \rangle \:\text{po rozwinięciu wykładnika sumy w iloczyn wykładników,} \\ & = \frac{1}{\sqrt{N}} \bigotimes_{k=1}^n \left(\vert0\rangle + e^{2 \pi i x /2^k } \vert1\rangle \right) \:\text{po uporządkowaniu sumy i iloczynów oraz rozwinięciu} \sum_{y=0}^{N-1} = \sum_{y_1=0}^{1}\sum_{y_2=0}^{1}\ldots\sum_{y_n=0}^{1} \\ & = \frac{1}{\sqrt{N}} \left(\vert0\rangle + e^{\frac{2\pi i}{2}x} \vert1\rangle\right) \otimes \left(\vert0\rangle + e^{\frac{2\pi i}{2^2}x} \vert1\rangle\right) \otimes \ldots \otimes \left(\vert0\rangle + e^{\frac{2\pi i}{2^{n-1}}x} \vert1\rangle\right) \otimes \left(\vert0\rangle + e^{\frac{2\pi i}{2^n}x} \vert1\rangle\right) \end{aligned}$$

2.4 Przykład: konstrukcja obwodu QFT dla 3 kubitów

Na podstawie powyższego opisu może nie być oczywiste, jak skonstruować obwód QFT. Na razie zauważmy po prostu, że oczekujemy, iż trzy kubity będą miały fazy ewoluujące z różnymi „szybkościami". Dokładne zrozumienie, jak przekłada się to na obwód, pozostawiamy czytelnikowi jako ćwiczenie. W pliku pdf z wykładem znajduje się wiele diagramów i przykładów. Dodatkowe materiały obejmują tę lekcję z kursu Fundamentals of quantum algorithms.

Będziemy demonstrować QFT wyłącznie przy użyciu symulatorów, w związku z tym nie będziemy korzystać ze struktury Qiskit patterns aż do momentu, gdy przejdziemy do kwantowej estymacji fazy.

# Prepare for 3 qubits circuit
qr = QuantumRegister(3)
cr = ClassicalRegister(3)
qc = QuantumCircuit(qr, cr)

qc.h(2)
qc.cp(pi / 2, 1, 2)  # Controlled-phase gate from qubit 1 to qubit 2
qc.cp(pi / 4, 0, 2)  # Controlled-phase gate from qubit 0 to qubit 2
qc.draw(output="mpl")

qc.h(1)
qc.cp(pi / 2, 0, 1)  # Controlled-phase gate from qubit 0 to qubit 1

qc.draw(output="mpl")

qc.h(0)

qc.draw(output="mpl")

qc.swap(0, 2)

qc.draw(output="mpl")

Jako przykład spróbujmy zastosować QFT do $|5\rangle$.

Najpierw potwierdzamy binarną reprezentację liczby całkowitej 5 i tworzymy obwód kodujący stan 5:

bin(5)

'0b101'

qc = QuantumCircuit(3)

qc.x(0)
qc.x(2)
qc.draw(output="mpl")

Sprawdzamy stany kwantowe przy użyciu symulatora Aer:

statevector = Statevector(qc)
plot_bloch_multivector(statevector)

Na koniec dodajemy QFT i oglądamy końcowy stan naszych kubitów:

qc.h(2)
qc.cp(pi / 2, 1, 2)
qc.cp(pi / 4, 0, 2)

qc.h(1)
qc.cp(pi / 2, 0, 1)

qc.h(0)

qc.swap(0, 2)

qc.draw(output="mpl")

statevector = Statevector(qc)
plot_bloch_multivector(statevector)

Widzimy, że nasza funkcja QFT zadziałała poprawnie. Kubit 0 został obrócony o $\tfrac{5}{8}$ pełnego obrotu, kubit 1 o $\tfrac{10}{8}$ pełnego obrotu (co odpowiada $\tfrac{1}{4}$ pełnego obrotu), a kubit 2 o $\tfrac{20}{8}$ pełnego obrotu (co odpowiada $\tfrac{1}{2}$ pełnego obrotu).

2.5 Ćwiczenie: QFT

(1) Zaimplementuj QFT dla 4 kubitów.

##your code goes here##

(2) Zastosuj QFT do $|14\rangle$, zasymuluj i narysuj wektor stanu na sferze Blocha.

##your code goes here##

Rozwiązanie ćwiczenia: QFT

(1)

qr = QuantumRegister(4)
cr = ClassicalRegister(4)
qc = QuantumCircuit(qr, cr)

qc.h(3)
qc.cp(pi / 2, 2, 3)
qc.cp(pi / 4, 1, 3)
qc.cp(pi / 8, 0, 3)

qc.h(2)
qc.cp(pi / 2, 1, 2)
qc.cp(pi / 4, 0, 2)

qc.h(1)
qc.cp(pi / 2, 0, 1)

qc.h(0)

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

qc.draw(output="mpl")

(2)

bin(14)

'0b1110'

qc = QuantumCircuit(4)

qc.x(1)
qc.x(2)
qc.x(3)
qc.draw("mpl")

qc.h(3)
qc.cp(pi / 2, 2, 3)
qc.cp(pi / 4, 1, 3)
qc.cp(pi / 8, 0, 3)

qc.h(2)
qc.cp(pi / 2, 1, 2)
qc.cp(pi / 4, 0, 2)

qc.h(1)
qc.cp(pi / 2, 0, 1)
qc.h(0)

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

qc.draw(output="mpl")

statevector = Statevector(qc)
plot_bloch_multivector(statevector)

3. Podstawy kwantowej estymacji fazy (QPE)

Mając operację unitarną $U$, QPE estymuje $\theta$ w $U\vert\psi \rangle =e^{\boldsymbol{2\pi i} \theta }|\psi \rangle$, ponieważ $U$ jest unitarna, wszystkie jej wartości własne mają normę 1.

3.1 Procedura

1. Konfiguracja

$\vert\psi\rangle$ znajduje się w jednym zestawie rejestrów kubitów. Dodatkowy zestaw $n$ kubitów tworzy rejestr liczący, w którym będziemy przechowywać wartość $2^n\theta$:

$$|\psi_0\rangle = \lvert 0 \rangle^{\otimes n} \lvert \psi \rangle$$

2. Superpozycja

Zastosuj $n$-bitową operację bramki Hadamard $H^{\otimes n}$ na rejestrze liczącym:

$$|\psi_1\rangle = {\frac {1}{2^{\frac {n}{2}}}}\left(|0\rangle +|1\rangle \right)^{\otimes n} \lvert \psi \rangle$$

3. Kontrolowane operacje unitarne

Musimy wprowadzić kontrolowaną operację unitarną $CU$, która stosuje operator unitarny $U$ na rejestrze docelowym tylko wtedy, gdy odpowiadający mu bit sterujący ma wartość $|1\rangle$. Ponieważ $U$ jest operatorem unitarnym z wektorem własnym $|\psi\rangle$ takim, że $U|\psi \rangle =e^{\boldsymbol{2\pi i} \theta }|\psi \rangle$, oznacza to:

$$U^{2^{j}}|\psi \rangle =U^{2^{j}-1}U|\psi \rangle =U^{2^{j}-1}e^{2\pi i\theta }|\psi \rangle =\cdots =e^{2\pi i2^{j}\theta }|\psi \rangle$$

3.2 Przykład: QPE dla T-gate

Użyjmy bramki $T$ jako przykładu QPE i oszacujmy jej fazę $\theta$.

$$T|1\rangle = \begin{pmatrix} 1 & 0\\ 0 & e^\frac{i\pi}{4}\\ \end{pmatrix} \begin{pmatrix} 0\\ 1\\ \end{pmatrix} = e^\frac{i\pi}{4}|1\rangle$$

Oczekujemy, że otrzymamy:

$$\theta = \frac{1}{8}$$

gdzie

$$T|1\rangle = e^{2i\pi\theta}|1\rangle$$

Inicjalizujemy $\vert\psi\rangle = \vert1\rangle$ jako wektor własny bramki $T$, stosując bramkę $X$:

qpe = QuantumCircuit(4, 3)
qpe.x(3)
qpe.draw(output="mpl")

Następnie stosujemy bramki Hadamard do kubitów liczących:

for qubit in range(3):
    qpe.h(qubit)
qpe.draw(output="mpl")

Wykonujemy kontrolowane operacje unitarne:

repetitions = 1
for counting_qubit in range(3):
    for i in range(repetitions):
        qpe.cp(pi / 4, counting_qubit, 3)  # This is C-U
    repetitions *= 2
qpe.draw(output="mpl")

Stosujemy odwrotną kwantową transformatę Fouriera, aby przekształcić stan rejestru liczącego, a następnie mierzymy rejestr liczący:

from qiskit.circuit.library import QFT

# Apply inverse QFT
qpe.append(QFT(3, inverse=True), [0, 1, 2])
qpe.draw(output="mpl")

for n in range(3):
    qpe.measure(n, n)
qpe.draw(output="mpl")

Możemy przeprowadzić symulację, korzystając z symulatora Aer:

aer_sim = AerSimulator()
shots = 2048

pm = generate_preset_pass_manager(backend=aer_sim, optimization_level=1)
t_qpe = pm.run(qpe)

sampler = Sampler(mode=aer_sim)
job = sampler.run([t_qpe], shots=shots)
result = job.result()
answer = result[0].data.c.get_counts()

plot_histogram(answer)

Widzimy, że otrzymujemy jeden wynik (001) z pewnością, co w systemie dziesiętnym odpowiada: 1. Teraz musimy podzielić nasz wynik (1) przez $2^n$, aby otrzymać $\theta$:

$$\theta = \frac{1}{2^3} = \frac{1}{8}$$

Algorytm Shora pozwala nam na faktoryzację liczby poprzez zbudowanie obwodu z nieznanym $\theta$ i uzyskanie $\theta$.

3.3 Ćwiczenie

Oszacuj $\theta = 1/3$, używając 3 kubitów do zliczania oraz jednego kubitu dla wektora własnego.

$P|1\rangle = e^{i\lambda}|1\rangle$. Ponieważ chcemy zaimplementować $U|1\rangle = e^{2\pi i \tfrac{1}{3}}|1\rangle$, musimy ustawić $\lambda = \tfrac{2 \pi}{3}$.

##your code goes here##

Rozwiązanie ćwiczenia: $\theta = 1/3$

# Create and set up circuit
qpe = QuantumCircuit(4, 3)

# Apply H-Gates to counting qubits:
for qubit in range(3):
    qpe.h(qubit)

# Prepare our eigenstate |psi>:
qpe.x(3)

# Do the controlled-U operations:
angle = 2 * pi / 3
repetitions = 1
for counting_qubit in range(3):
    for i in range(repetitions):
        qpe.cp(angle, counting_qubit, 3)
    repetitions *= 2

# Do the inverse QFT:
qpe.append(QFT(3, inverse=True), [0, 1, 2])

for n in range(3):
    qpe.measure(n, n)

qpe.draw(output="mpl")

aer_sim = AerSimulator()
shots = 4096

pm = generate_preset_pass_manager(backend=aer_sim, optimization_level=1)
t_qpe = pm.run(qpe)

sampler = Sampler(mode=aer_sim)
job = sampler.run([t_qpe], shots=shots)
result = job.result()
answer = result[0].data.c.get_counts()

plot_histogram(answer)

4. Wykonanie z użyciem prymitywu Sampler z Qiskit Runtime

Wykonamy QPE na rzeczywistym urządzeniu kwantowym, przechodząc przez 4 kroki wzorców Qiskit.

1. Mapowanie problemu na format natywny dla komputera kwantowego
2. Optymalizacja obwodów
3. Wykonanie docelowego obwodu
4. Post-processing wyników

from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import Sampler

# Loading your IBM Quantum® account(s)

service = QiskitRuntimeService()

4.1 Krok 1: Mapowanie problemu na obwody i operatory kwantowe

qpe = QuantumCircuit(4, 3)
qpe.x(3)
for qubit in range(3):
    qpe.h(qubit)
repetitions = 1
for counting_qubit in range(3):
    for i in range(repetitions):
        qpe.cp(pi / 4, counting_qubit, 3)  # This is C-U
    repetitions *= 2
qpe.append(QFT(3, inverse=True), [0, 1, 2])
for n in range(3):
    qpe.measure(n, n)
qpe.draw(output="mpl")

backend = service.least_busy(simulator=False, operational=True, min_num_qubits=4)

print(backend)

<IBMBackend('ibm_strasbourg')>

4.2 Krok 2: Optymalizacja pod docelowy sprzęt

# Transpile the circuit into basis gates executable on the hardware
pm = generate_preset_pass_manager(backend=backend, optimization_level=2)
qc_compiled = pm.run(qpe)

qc_compiled.draw("mpl", idle_wires=False)

4.3 Krok 3: Wykonanie na docelowym sprzęcie

real_sampler = Sampler(mode=backend)
job = real_sampler.run([qc_compiled], shots=1024)
job_id = job.job_id()
print("job id:", job_id)

job id: d13p4zb5z6q00087ag00

job = service.job(job_id)  # Input your job-id between the quotations
job.status()

'DONE'

result_real = job.result()
print(result_real)

PrimitiveResult([SamplerPubResult(data=DataBin(c=BitArray(<shape=(), num_shots=1024, num_bits=3>)), metadata={'circuit_metadata': {}})], metadata={'execution': {'execution_spans': ExecutionSpans([DoubleSliceSpan(<start='2025-06-09 22:39:00', stop='2025-06-09 22:39:00', size=1024>)])}, 'version': 2})

4.4 Krok 4: Przetwarzanie końcowe wyników

from qiskit.visualization import plot_histogram

plot_histogram(result_real[0].data.c.get_counts())

# See the version of Qiskit
import qiskit

qiskit.__version__

'2.0.2'