Simule sistemas quânticos com ruído usando propagação de Pauli

Neste guia, usamos o pacote pauli-prop para simular classicamente a dinâmica temporal de um modelo de Ising de campo transversal (TFIM) ruidoso de nove qubits em uma rede quadrada 3x3. Usamos instruções PauliLindbladError para definir um canal de ruído, $\Lambda$, que atua sobre um conjunto de camadas de emaranhamento, $\mathcal{U}$. Em seguida, propagamos o observável, $O$, de trás para frente através do circuito ruidoso e estimamos valores esperados para uma variedade de modelos de ruído, bem como para o caso sem ruído.

À medida que o observável é propagado para trás através do circuito, cada canal de ruído, $\Lambda_k$, associado à camada de emaranhamento, $\mathcal{U}_k$, amortece os termos de Pauli em $O$ que anticomutam com seus geradores de Pauli-Lindblad. Especificamente, se $G_{k,i}$ é um gerador de Pauli de $\Lambda_k$ com taxa $\gamma_{k,i}$, então um termo de Pauli, $P$, em $O$ se transforma como: $c_P \mapsto c_P e^{-2\gamma_{k,i}} \quad \text{if } \{P, G_{k,i}\}=0$, onde $c_P$ é o coeficiente de $P$. Quando $O$ tiver sido propagado até o início do circuito, o valor esperado em relação ao estado zero, $|0\rangle^{\otimes N}$, pode ser calculado trivialmente somando os coeficientes de cada termo diagonal em $O$ (termos que contêm $Z$ ou $I$ em todos os qubits).

Fluxo de trabalho:

• Especifique a rede do TFIM e use coloração de arestas para identificar um conjunto mínimo de camadas de emaranhamento
• Gere modelos de ruído sintéticos, $\Lambda_k$, para cada camada de emaranhamento única, $U_k$
  • Crie modelos de ruído de várias escalas para estudar o impacto do ruído das portas no sistema
• Crie circuitos quânticos sem ruído e com ruído para as várias profundidades e escalas de ruído de interesse
  • Nos circuitos com ruído, instruções PauliLindbladError são inseridas antes de cada camada de emaranhamento
• Use a propagação de Pauli para simular valores esperados exatos do sistema em várias profundidades
  • Para nove qubits, isso é feito deixando $O$ crescer até $4^9$ termos, cobrindo todo o espaço de Pauli
• Use a propagação de Pauli para simular valores esperados com ruído
• Observe como o aumento do ruído das portas degrada a precisão do modelo quântico

Gere uma rede quadrada 3x3 e encontre uma coloração com 4 cores nas arestas

Os vértices do grafo representam qubits, e as arestas representam uma conexão entre dois qubits. A coloração das arestas corresponde às camadas de emaranhamento únicas no circuito quântico, de modo que portas em conexões associadas a cores diferentes não podem ser aplicadas simultaneamente.

Identificar um conjunto mínimo de camadas de emaranhamento únicas costuma ser importante para implementar protocolos eficientes de aprendizado de ruído, pois o ruído de cada camada deve ser aprendido de forma independente. Quanto mais camadas precisamos aprender, mais shots precisamos coletar da QPU. Nesta demonstração, usamos as informações das camadas para construir circuitos com ruído e injetar instruções PauliLindbladError de qiskit-aer antes de cada camada de emaranhamento, para modelar o ruído das portas da QPU.

# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util

_needed = {"matplotlib": "matplotlib", "numpy": "numpy", "pauli_prop": "pauli-prop", "qiskit": "qiskit", "qiskit_addon_utils": "qiskit-addon-utils", "qiskit_aer": "qiskit-aer", "rustworkx": "rustworkx"}
_missing = [pip for module, pip in _needed.items()
            if importlib.util.find_spec(module) is None]
# One at a time, so a package that fails to install does not block the others
for _pip in _missing:
    %pip install -q {_pip}
if not _missing:
    print("\u2713 All required packages are installed")

from collections import defaultdict

import numpy as np
from qiskit.transpiler import CouplingMap
from qiskit_addon_utils.coloring import auto_color_edges

# Define rectangular square-lattice on 20 qubits
num_rows = 3
num_cols = 3
num_qubits = num_rows * num_cols

coupling_map = CouplingMap.from_grid(
    num_rows=num_rows, num_columns=num_cols, bidirectional=False
)

# Create mapping from color to edge list
coloring = auto_color_edges(coupling_map.get_edges())
color_to_edge = defaultdict(list)
for edge, color in coloring.items():
    color_to_edge[color].append(edge)

from rustworkx import PyDiGraph
from rustworkx.visualization import graphviz_draw

# Inspect graph coupling and unique entangling layers
print(
    f"The circuit will have {num_qubits} qubits and {len(color_to_edge)} unique entangling layers."
)
sq_lattice = PyDiGraph()
sq_lattice.extend_from_weighted_edge_list(
    [
        (source, target, color)
        for ((source, target), color) in coloring.items()
    ]
)

def color_edge_4color(edge):
    color_dict = {0: "red", 1: "green", 2: "blue", 3: "orange"}
    return {"color": color_dict[edge]}

graphviz_draw(sq_lattice, edge_attr_fn=color_edge_4color, method="neato")

The circuit will have 9 qubits and 4 unique entangling layers.

Gere modelos de ruído sintéticos

Antes de criar os circuitos quânticos, geramos um modelo de ruído (instância de PauliLindbladError) para cada uma das camadas de emaranhamento. Vamos incorporá-los como instruções em nossos circuitos quânticos mais adiante. Para cada camada, geramos canais de ruído de várias escalas. Especificamente, geramos modelos de ruído com Error Per Layered Gate (EPLG) de aproximadamente .0004, .0008, .0012, .0016, e .002.

from qiskit.quantum_info import SparsePauliOp, pauli_basis
from qiskit_aer.noise import PauliLindbladError

# Pauli-Lindblad noise parameters
seed = 1764
target_EPLGs = [0.0004, 0.0008, 0.0012, 0.0016, 0.002]

def generate_random_pauli_lindblad_noise(
    edges,
    num_qubits: int | None = None,
    noise_scale: float = 1e-3,
    seed: int | None = None,
) -> PauliLindbladError:
    """Generate random Pauli-Lindblad noise over the full Pauli basis."""
    if num_qubits is None:
        num_qubits = np.max(edges)

    basis_paulis = [p for p in pauli_basis(2) if np.sum(p.x + p.z)]
    basis_paulis = SparsePauliOp.from_sparse_list(
        [
            (pauli.to_label(), edge, 1)
            for pauli in basis_paulis
            for edge in edges
        ],
        num_qubits=num_qubits,
    )
    basis_paulis = basis_paulis.simplify()
    basis_paulis = basis_paulis.paulis

    rng = np.random.default_rng(seed=seed)
    rates = rng.random(len(basis_paulis)) * noise_scale

    return PauliLindbladError(generators=basis_paulis, rates=rates)

num_generators = (
    (num_rows * num_cols)
    + (num_rows - 1) * num_cols
    + num_rows * (num_cols - 1)
)
noise_scales = [
    EPLG * (num_rows * num_cols) / num_generators for EPLG in target_EPLGs
]
noise_models_per_EPLG = [
    [
        generate_random_pauli_lindblad_noise(
            color_to_edge[color],
            num_qubits=num_qubits,
            noise_scale=noise_scale,
            seed=seed,
        )
        for color in range(len(color_to_edge))
    ]
    for noise_scale in noise_scales
]

Crie os circuitos quânticos

Nesta demonstração, simulamos a dinâmica temporal de um modelo de Ising de campo transversal (TFIM) para números crescentes de passos de Trotter (1 a 10 passos). Para cada uma das 10 profundidades de circuito, simulamos o efeito do ruído das portas, dados modelos de ruído de várias escalas (EPLGs = .0004, .0008, .0012, .0016, .002). O ruído é inserido no QuantumCircuit como uma instrução PauliLindbladError do Qiskit Aer. O Hamiltoniano considerado é:

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

onde $J>0$ descreve o acoplamento de spins vizinhos mais próximos, $i<j$, e $h$ é o campo transversal global.

Aqui implementamos o Hamiltoniano evoluído no tempo em várias escalas de tempo e de ruído. Criamos 60 circuitos no total: 10 circuitos sem ruído que variam na profundidade de Trotter e 50 circuitos com ruído para as 10 profundidades de Trotter em cinco escalas de ruído. Dado um grafo de conectividade, o modelo é parametrizado por algumas variáveis:

• num_steps: O número de passos de Trotter
• J: Intensidade do acoplamento entre sítios conectados
• h: Intensidade do campo magnético externo
• dt: Variação do tempo ao longo de um passo de Trotter
• initial_state_angle: Uma excitação inicial, $R_y(\theta)$, a ser aplicada uniformemente a todos os qubits

from typing import Any

from qiskit import QuantumCircuit

# Ising model parameters
num_steps = 10
J = -1.0
dt = 0.25 / abs(J)
h = 2.0 * abs(J)
initial_state_angle = np.pi / 18.0
rx_angle = 2.0 * h * dt
rzz_angle = 2.0 * J * dt

def generate_ising_circuit(
    num_qubits: int,
    num_steps: int,
    rx_angle: float,
    rzz_angle: float,
    coloring: dict[Any, list[tuple[int, int]]],
    layer_noise_models: list[PauliLindbladError] | None = None,
    initial_state_angle: float | None = None,
) -> QuantumCircuit:
    """Generate a quantum circuit implementing a transverse-field Ising model"""
    qc = QuantumCircuit(num_qubits)
    if initial_state_angle:
        qc.ry(initial_state_angle, range(num_qubits))
    qc.rx(rx_angle / 2, range(num_qubits))
    for i in range(num_steps):
        for j, layer in enumerate(coloring):
            edges = coloring[layer]
            if layer_noise_models:
                qc.append(layer_noise_models[j], qargs=range(num_qubits))
            for edge in edges:
                qc.rzz(rzz_angle, *edge)
        if i == num_steps - 1:
            qc.rx(rx_angle / 2, range(num_qubits))
        else:
            qc.rx(rx_angle, range(num_qubits))
    return qc

# Create the noiseless and noisy circuits
noiseless_circs = []
noisy_circs = []
for steps in range(1, num_steps + 1):
    noiseless_circs.append(
        generate_ising_circuit(
            num_qubits,
            steps,
            rx_angle,
            rzz_angle,
            color_to_edge,
            initial_state_angle=initial_state_angle,
        )
    )
    noisy_circs_per_step = []
    for noise_models in noise_models_per_EPLG:
        noisy_circs_per_step.append(
            generate_ising_circuit(
                num_qubits,
                steps,
                rx_angle,
                rzz_angle,
                color_to_edge,
                layer_noise_models=noise_models,
                initial_state_angle=initial_state_angle,
            )
        )
    noisy_circs.append(noisy_circs_per_step)
print(
    f"{num_steps} noiseless and {num_steps * len(target_EPLGs)} noisy Trotter circuits generated. {num_steps} different depths across {len(target_EPLGs)} different noise models"
)
print("\nBelow: Initial state and one noisy Trotter step.")
noisy_circs[0][0].draw("mpl", fold=-1)

10 noiseless and 50 noisy Trotter circuits generated. 10 different depths across 5 different noise models

Below: Initial state and one noisy Trotter step.

Especifique o observável e execute as simulações

Nesta demonstração, simulamos valores esperados do correlador médio de dois sítios:

$\langle O \rangle = \langle Z_{tot}^2(s) \rangle = \frac{1}{N^2}\sum \langle \Psi(\theta)|(\mathscr{U}^{\dagger})^sZ_jZ_k(\mathscr{U})^s|\Psi(\theta) \rangle$

onde $\Psi(\theta)$ corresponde a uma rotação uniforme $R_y(\theta)$ em todos os qubits, $\mathscr{U}^s$ descreve $s$ camadas de Trotter, e $(j,k)$ indexa todos os pares de vértices conectados da rede.

Por fim, usamos pauli_prop para simular os valores esperados do observável para cada um dos circuitos. Nesta demonstração de nove qubits, realizamos todas as simulações de forma exata. Nenhum truncamento da propagação de Pauli será realizado; portanto, as diferenças nos valores esperados entre os diferentes modelos de ruído podem ser atribuídas inteiramente ao erro das portas. O processo de simulação é feito em quatro etapas:

• Evolua as portas Clifford do circuito para o início do circuito usando pauli_prop.evolve_through_cliffords
• Propague o observável através da parte não Clifford do circuito usando pauli_prop.propagate_through_circuit
  • Realizamos simulações exatas permitindo que o observável cresça até o tamanho de todo o espaço de Pauli, $4^9$
• Propague o observável evoluído através da parte Clifford do circuito usando SparsePauliOp.evolve do Qiskit
• Estime o valor esperado em relação ao estado zero, $|0\rangle^{\otimes N}$, somando os coeficientes de cada termo diagonal em $O$ (termos que contêm $Z$ ou $I$ em todos os qubits)

import time

from pauli_prop import evolve_through_cliffords, propagate_through_circuit
from qiskit.quantum_info import Pauli

# Average ZZ-correlator observable
id_pauli = Pauli("I" * num_qubits)
observable = 2 * SparsePauliOp(
    [
        id_pauli.dot(Pauli("ZZ"), [i, j])
        for i in range(num_qubits)
        for j in range(i + 1, num_qubits)
    ]
)
observable /= num_qubits**2

# Pauli propagation parameters
max_terms = 4**num_qubits  # Exact propagation
atol = 1e-12

# Run simulations
exact_evs = []
noisy_evs = [[] for _ in range(len(target_EPLGs))]
st = time.perf_counter()
for i, noiseless_circ in enumerate(noiseless_circs):
    cliff, non_cliff = evolve_through_cliffords(noiseless_circ)
    evolved_obs = propagate_through_circuit(
        observable, non_cliff, max_terms=max_terms, atol=atol, frame="h"
    )[0]
    evolved_obs.paulis = evolved_obs.paulis.evolve(cliff, frame="h")
    exact_evs.append(
        float(evolved_obs.coeffs[~evolved_obs.paulis.x.any(axis=1)].sum())
    )
    for j in range(len(target_EPLGs)):
        noisy_circ = noisy_circs[i][j]
        cliff, non_cliff = evolve_through_cliffords(noisy_circ)
        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")
        noisy_evs[j].append(
            float(evolved_obs.coeffs[~evolved_obs.paulis.x.any(axis=1)].sum())
        )
print(
    f"Ran {len(noiseless_circs)} noiseless and {len(target_EPLGs) * num_steps} noisy simulations in {int(time.perf_counter() - st)}s."
)

Ran 10 noiseless and 50 noisy simulations in 103s.

Observe o efeito do erro das portas no modelo

Lembre-se de que, como este é um experimento de nove qubits, a rotina de propagação de Pauli é exata, e todo o erro nos gráficos de ruído pode ser atribuído ao erro das portas.

import matplotlib.pyplot as plt

xs = range(1, num_steps + 1)
plt.plot(xs, exact_evs, label="Noiseless", color="black", marker="o")
colors = [".3", ".4", ".5", ".6", ".7"]
for i, evs in enumerate(noisy_evs):
    plt.plot(
        xs,
        evs,
        label=f"{target_EPLGs[i]} EPLG",
        linestyle="--",
        color=colors[i],
        marker="o",
    )
plt.xlabel("# Trotter steps")
plt.ylabel(r"$\langle Z_{tot}^2 \rangle$")
plt.legend()
plt.show()