title: "Avvio rapido" description: "Guida rapida al pacchetto di componenti aggiuntivi pauli-prop per Qiskit"

Avvio rapido

In questa guida utilizziamo il pacchetto pauli-prop per simulare in modo classico le dinamiche temporali di un modello di Ising “kicked” a 10 qubit su una catena di spin di tipo “ 1D ”.

Preparare i dati di input per la propagazione di Pauli

L'Hamiltoniano considerato è:

$H = -J\sum\limits_{\langle i,j \rangle} Z_iZ_j + h\sum\limits_iX_i$

dove $J>0$ descrive l'accoppiamento tra gli spin più vicini, $i<j$, e $h$ è il campo trasversale globale. Una decomposizione di Trotter di primo ordine dell’operatore evoluto nel tempo verrà implementata come circuito quantistico, $U$, in $20$ passi di Trotter. La costante di accoppiamento, $J$, sarà fissata a $J=-\frac{\pi}{2}$, mentre $h$ sarà fissata a $\frac{\pi}{6}$. Le interazioni $ZZ$ saranno implementate utilizzando i gate di Clifford ( $CX$, $Sdg$, $\sqrt{Y}$ ).

Implementiamo l’evoluzione temporale “trotterizzata” come circuito quantistico e utilizziamo l’ $\frac{\pi}{6}$ e per le rotazioni non-Clifford attorno all’asse x. Quanto più questi angoli si discostano dagli angoli di Clifford (ad esempio, $\theta=n\frac{\pi}{2}, n \in \mathbb{Z}$ ), tanto più difficile sarà simulare il sistema con i metodi di propagazione di Pauli.

Per la scelta dell’osservabile, consideriamo la magnetizzazione media di un singolo sito, $\frac{1}{N} \sum_{i=1}^{N} \langle z_i \rangle$, dove $N$ è il numero di spin.

import numpy as np
from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit.transpiler import CouplingMap

num_qubits = 10
coupling_map = CouplingMap.from_line(num_qubits, bidirectional=False)

# Num Trotter steps
num_steps = 20
theta_rx = np.pi / 6

# Average single-site magnetization
observable = (
    SparsePauliOp(
        [
            "I" * iq + "Z" + "I" * (num_qubits - iq - 1)
            for iq in range(num_qubits)
        ]
    )
    / num_qubits
)

# Create the Trotter circuit
num_qubits = 10
num_steps = 20
theta_rx = np.pi / 6
circuit = QuantumCircuit(num_qubits)
edges = CouplingMap.from_line(num_qubits, bidirectional=False).get_edges()
for _ in range(num_steps):
    circuit.rx(theta_rx, [i for i in range(num_qubits)])
    for edge in edges:
        circuit.sdg(edge)
        circuit.ry(np.pi / 2, edge[1])
        circuit.cx(edge[0], edge[1])
        circuit.ry(-np.pi / 2, edge[1])
circuit.draw("mpl", fold=-1)

Output:

Simulare l'evoluzione temporale del sistema con la propagazione di Pauli

Una volta definiti il circuito ( $U$ ) e l’osservabile ( $O$ ), possiamo simulare facilmente il sistema in pochi passaggi:

• Separare $U$ nelle sue parti relative a Clifford, $C$, e non relative a Clifford, $P$, in modo tale che $U=PC$ utilizzando evolve_through_cliffords
• Applicare la funzione " $O$ " alla funzione " $P$ ", ottenendo una nuova funzione, " $O^\prime$ ", utilizzando pauli_prop.propagate_through_circuit
• Evolvere l' $O^\prime$ e lungo la parte Clifford del circuito utilizzando il supporto per l'evoluzione Clifford integrato in Qiskit
• Approssimare il valore atteso come $\langle0|O^\prime|0\rangle \approx \langle0|U^\dagger OU|0\rangle$ sommando i coefficienti in $O^\prime$ associati ai termini di Pauli completamente diagonali (termini di Pauli contenenti o I Z su tutti i qubit). Ricordate che si tratta di un'approssimazione, poiché abbiamo troncato i termini provenienti dall' $O^\prime$ e mentre li propagavamo attraverso la parte non-Clifford del circuito.

import time

from pauli_prop import evolve_through_cliffords, propagate_through_circuit

cliff, non_cliff = evolve_through_cliffords(circuit)

max_terms_list = [10**i for i in range(8)]
approx_evs = []
durations = []
for max_terms in max_terms_list:
    st = time.perf_counter()
    evolved_obs = propagate_through_circuit(
        observable, non_cliff, max_terms=max_terms, atol=1e-12, frame="h"
    )[0]
    evolved_obs.paulis = evolved_obs.paulis.evolve(cliff, frame="h")
    durations.append(time.perf_counter() - st)
    approx_evs.append(
        float(evolved_obs.coeffs[~evolved_obs.paulis.x.any(axis=1)].sum())
    )

Man mano che si eseguono calcoli più complessi, le approssimazioni del valore atteso diventano più precise. In questo esempio, saturiamo l'intero spazio di Pauli a circa $4^{10}\approx10^6$, il che si riflette nell'appiattimento della curva tra gli ultimi due punti.

Sebbene il grafico riportato di seguito mostri una convergenza monotona, le simulazioni di propagazione di Pauli generalmente non convergono in modo monotono. Non è raro osservare un andamento “irregolare” in questo tipo di grafici.

import matplotlib.pyplot as plt
from qiskit_aer import AerSimulator

sim_circ = circuit.copy()
sim_circ.save_statevector()
backend = AerSimulator(method="statevector")
psi = backend.run(sim_circ).result().data()["statevector"]
exact_ev = psi.expectation_value(observable)

ax1 = plt.gca()
ax1.plot(max_terms_list, approx_evs, marker="o", label="Approximate")
ax1.axhline(exact_ev, linestyle="--", color="green", label="Exact")
ax1.set_xscale("log")
ax1.set_xlabel("# terms kept")
ax1.set_ylabel(r"$\frac{1}{N} \sum_{i=1}^{N} \langle z_i \rangle$")

ax2 = ax1.twinx()
ax2.plot(
    max_terms_list, durations, marker=".", label="Runtime", color="orange"
)
ax2.set_ylabel("Runtime (s)", color="orange")
ax2.set_yscale("log")

handles1, labels1 = ax1.get_legend_handles_labels()
handles2, labels2 = ax2.get_legend_handles_labels()
ax1.legend(handles1 + handles2, labels1 + labels2, loc="lower right")

plt.title(f"Simulating {num_steps}-step 1D Ising Model")

Output:

Text(0.5, 1.0, 'Simulating 20-step 1D Ising Model')