Modelo de Ising con campo transversal con la gestión de rendimiento de Q-CTRL

Estimación de uso: 2 minutos en un procesador Heron r2. (NOTA: Esto es solo una estimación. Su tiempo de ejecución puede variar.)

Contexto

El modelo de Ising con campo transversal (TFIM, por sus siglas en inglés) es importante para estudiar el magnetismo cuántico y las transiciones de fase. Describe un conjunto de espines dispuestos en una red, donde cada espín interactúa con sus vecinos mientras también es influenciado por un campo magnético externo que impulsa las fluctuaciones cuánticas.

Un enfoque común para simular este modelo es utilizar la descomposición de Trotter para aproximar el operador de evolución temporal, construyendo circuitos que alternan entre rotaciones de un solo cúbit y acoplamientos de dos cúbits. Sin embargo, esta simulación en hardware real es desafiante debido al ruido y la decoherencia, lo que conduce a desviaciones de la dinámica verdadera. Para superar esto, utilizamos las herramientas de supresión de errores y gestión de rendimiento Fire Opal de Q-CTRL, ofrecidas como una función de Qiskit (consulta la documentación de Fire Opal). Fire Opal optimiza automáticamente la ejecución de circuitos aplicando desacoplamiento dinámico, diseño avanzado de layout, enrutamiento y otras técnicas de supresión de errores, todo con el objetivo de reducir el ruido. Con estas mejoras, los resultados del hardware se alinean más estrechamente con las simulaciones sin ruido, y así podemos estudiar la dinámica de magnetización del TFIM con mayor fidelidad.

En este tutorial:

• Construiremos el hamiltoniano del TFIM en un grafo de triángulos de espines conectados
• Simularemos la evolución temporal con circuitos trotterizados a diferentes profundidades
• Calcularemos y visualizaremos las magnetizaciones de un solo cúbit $\langle Z_i \rangle$ a lo largo del tiempo
• Compararemos simulaciones de referencia con resultados de ejecuciones en hardware utilizando la gestión de rendimiento Fire Opal de Q-CTRL

Descripción general

El modelo de Ising con campo transversal (TFIM) es un modelo de espín cuántico que captura las características esenciales de las transiciones de fase cuánticas. El hamiltoniano se define como:

$$H = -J \sum_{i} Z_i Z_{i+1} - h \sum_{i} X_i$$

donde $Z_i$ y $X_i$ son operadores de Pauli que actúan sobre el cúbit $i$, $J$ es la intensidad del acoplamiento entre espines vecinos, y $h$ es la intensidad del campo magnético transversal. El primer término representa las interacciones ferromagnéticas clásicas, mientras que el segundo introduce fluctuaciones cuánticas a través del campo transversal. Para simular la dinámica del TFIM, se utiliza una descomposición de Trotter del operador de evolución unitaria $e^{-iHt}$, implementada a través de capas de compuertas RX y RZZ basadas en un grafo personalizado de triángulos de espines conectados. La simulación explora cómo la magnetización $\langle Z \rangle$ evoluciona con el aumento de pasos de Trotter.

El rendimiento de la implementación propuesta del TFIM se evalúa comparando simulaciones sin ruido con backends ruidosos. Las características de ejecución mejorada y supresión de errores de Fire Opal se utilizan para mitigar el efecto del ruido en hardware real, proporcionando estimaciones más confiables de los observables de espín como $\langle Z_i \rangle$ y los correladores $\langle Z_i Z_j \rangle$.

Requisitos

Antes de comenzar este tutorial, asegúrate de tener instalado lo siguiente:

• Qiskit SDK v1.4 o posterior, con soporte de visualización
• Qiskit Runtime v0.40 o posterior (pip install qiskit-ibm-runtime)
• Qiskit Functions Catalog v0.9.0 (pip install qiskit-ibm-catalog)
• Fire Opal SDK v9.0.2 o posterior (pip install fire-opal)
• Q-CTRL Visualizer v8.0.2 o posterior (pip install qctrl-visualizer)

Configuración

Primero, autentíquese utilizando su clave API de IBM Quantum. Luego, seleccione la función de Qiskit de la siguiente manera. (Este código asume que ya ha guardado su cuenta en su entorno local.)

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

from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit import QuantumCircuit
from qiskit_ibm_catalog import QiskitFunctionsCatalog
from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import SamplerV2 as Sampler
from qiskit.quantum_info import SparsePauliOp
from qiskit_aer import AerSimulator

import numpy as np
import networkx as nx
import matplotlib.pyplot as plt
import qctrlvisualizer as qv

catalog = QiskitFunctionsCatalog(channel="ibm_quantum_platform")

# Access Function
perf_mgmt = catalog.load("q-ctrl/performance-management")

Paso 1: Asignar entradas clásicas a un problema cuántico

Generar el grafo del TFIM

Comenzamos definiendo la red de espines y los acoplamientos entre ellos. En este tutorial, la red se construye a partir de triángulos conectados dispuestos en una cadena lineal. Cada triángulo consta de tres nodos conectados en un lazo cerrado, y la cadena se forma vinculando un nodo de cada triángulo con el triángulo anterior.

La función auxiliar connected_triangles_adj_matrix construye la matriz de adyacencia para esta estructura. Para una cadena de $n$ triángulos, el grafo resultante contiene $2n+1$ nodos.

def connected_triangles_adj_matrix(n):
    """
    Generate the adjacency matrix for 'n' connected triangles in a chain.
    """
    num_nodes = 2 * n + 1
    adj_matrix = np.zeros((num_nodes, num_nodes), dtype=int)

    for i in range(n):
        a, b, c = i * 2, i * 2 + 1, i * 2 + 2  # Nodes of the current triangle

        # Connect the three nodes in a triangle
        adj_matrix[a, b] = adj_matrix[b, a] = 1
        adj_matrix[b, c] = adj_matrix[c, b] = 1
        adj_matrix[a, c] = adj_matrix[c, a] = 1

        # If not the first triangle, connect to the previous triangle
        if i > 0:
            adj_matrix[a, a - 1] = adj_matrix[a - 1, a] = 1

    return adj_matrix

Para visualizar la red que acabamos de definir, podemos graficar la cadena de triángulos conectados y etiquetar cada nodo. La función a continuación construye el grafo para un número elegido de triángulos y lo muestra.

def plot_triangle_chain(n, side=1.0):
    """
    Plot a horizontal chain of n equilateral triangles.
    Baseline: even nodes (0,2,4,...,2n) on y=0
    Apexes: odd nodes (1,3,5,...,2n-1) above the midpoint.
    """
    # Build graph
    A = connected_triangles_adj_matrix(n)
    G = nx.from_numpy_array(A)

    h = np.sqrt(3) / 2 * side
    pos = {}

    # Place baseline nodes
    for k in range(n + 1):
        pos[2 * k] = (k * side, 0.0)

    # Place apex nodes
    for k in range(n):
        x_left = pos[2 * k][0]
        x_right = pos[2 * k + 2][0]
        pos[2 * k + 1] = ((x_left + x_right) / 2, h)

    # Draw
    fig, ax = plt.subplots(figsize=(1.5 * n, 2.5))
    nx.draw(
        G,
        pos,
        ax=ax,
        with_labels=True,
        font_size=10,
        font_color="white",
        node_size=600,
        node_color=qv.QCTRL_STYLE_COLORS[0],
        edge_color="black",
        width=2,
    )
    ax.set_aspect("equal")
    ax.margins(0.2)
    plt.show()

    return G, pos

Para este tutorial utilizaremos una cadena de 20 triángulos.

n_triangles = 20
n_qubits = 2 * n_triangles + 1
plot_triangle_chain(n_triangles, side=1.0)
plt.show()

Colorear las aristas del grafo

Para implementar el acoplamiento espín-espín, es útil agrupar las aristas que no se superponen. Esto nos permite aplicar compuertas de dos cúbits en paralelo. Podemos hacer esto con un procedimiento simple de coloración de aristas [1], que asigna un color a cada arista de modo que las aristas que se encuentran en el mismo nodo se coloquen en grupos diferentes.

def edge_coloring(graph):
    """
    Takes a NetworkX graph and returns a list of lists where each inner list contains
    the edges assigned the same color.
    """
    line_graph = nx.line_graph(graph)
    edge_colors = nx.coloring.greedy_color(line_graph)

    color_groups = {}
    for edge, color in edge_colors.items():
        if color not in color_groups:
            color_groups[color] = []
        color_groups[color].append(edge)

    return list(color_groups.values())

Paso 2: Optimizar el problema para la ejecución en hardware cuántico

Generar circuitos trotterizados en grafos de espines

Para simular la dinámica del TFIM, construimos circuitos que aproximan el operador de evolución temporal.

$$U(t) = e^{-i H t}, \quad \text{where} \quad H = -J \sum_{\langle i,j \rangle} Z_i Z_j - h \sum_i X_i .$$

Utilizamos una descomposición de Trotter de segundo orden:

$$e^{-i H \Delta t} \approx e^{-i H_X \Delta t / 2}\, e^{-i H_Z \Delta t}\, e^{-i H_X \Delta t / 2},$$

donde $H_X = -h \sum_i X_i$ y $H_Z = -J \sum_{\langle i,j \rangle} Z_i Z_j$.

• El término $H_X$ se implementa con capas de rotaciones RX.
• El término $H_Z$ se implementa con capas de compuertas RZZ a lo largo de las aristas del grafo de interacción.

Los ángulos de estas compuertas están determinados por el campo transversal $h$, la constante de acoplamiento $J$ y el paso temporal $\Delta t$. Al apilar múltiples pasos de Trotter, generamos circuitos de profundidad creciente que aproximan la dinámica del sistema. Las funciones generate_tfim_circ_custom_graph y trotter_circuits construyen un circuito cuántico trotterizado a partir de un grafo de interacción de espines arbitrario.

def generate_tfim_circ_custom_graph(
    steps, h, J, dt, psi0, graph: nx.graph.Graph, meas_basis="Z", mirror=False
):
    """
    Generate a second order trotter of the form e^(a+b) ~ e^(b/2) e^a e^(b/2) for simulating a transverse field ising model:
    e^{-i H t} where the Hamiltonian H = -J \\sum_i Z_i Z_{i+1} + h \\sum_i X_i.

    steps: Number of trotter steps
    theta_x: Angle for layer of X rotations
    theta_zz: Angle for layer of ZZ rotations
    theta_x: Angle for second layer of X rotations
    J: Coupling between nearest neighbor spins
    h: The transverse magnetic field strength
    dt: t/total_steps
    psi0: initial state (assumed to be prepared in the computational basis).
    meas_basis: basis to measure all correlators in

    This is a second order trotter of the form e^(a+b) ~ e^(b/2) e^a e^(b/2)
    """
    theta_x = h * dt
    theta_zz = -2 * J * dt
    nq = graph.number_of_nodes()
    color_edges = edge_coloring(graph)
    circ = QuantumCircuit(nq, nq)
    # Initial state, for typical cases in the computational basis
    for i, b in enumerate(psi0):
        if b == "1":
            circ.x(i)
    # Trotter steps
    for step in range(steps):
        for i in range(nq):
            circ.rx(theta_x, i)
        if mirror:
            color_edges = [sublist[::-1] for sublist in color_edges[::-1]]
        for edge_list in color_edges:
            for edge in edge_list:
                circ.rzz(theta_zz, edge[0], edge[1])
        for i in range(nq):
            circ.rx(theta_x, i)

    # some typically used basis rotations
    if meas_basis == "X":
        for b in range(nq):
            circ.h(b)
    elif meas_basis == "Y":
        for b in range(nq):
            circ.sdg(b)
            circ.h(b)

    for i in range(nq):
        circ.measure(i, i)

    return circ

def trotter_circuits(G, d_ind_tot, J, h, dt, meas_basis, mirror=True):
    """
    Generates a sequence of Trotterized circuits, each with increasing depth.
    Given a spin interaction graph and Hamiltonian parameters, it constructs
    a list of circuits with 1 to d_ind_tot Trotter steps

    G: Graph defining spin interactions (edges = ZZ couplings)
    d_ind_tot: Number of Trotter steps (maximum depth)
    J: Coupling between nearest neighboring spins
    h: Transverse magnetic field strength
    dt: (t / total_steps
    meas_basis: Basis to measure all correlators in
    mirror: If True, mirror the Trotter layers
    """
    qubit_count = len(G)
    circuits = []
    psi0 = "0" * qubit_count

    for steps in range(1, d_ind_tot + 1):
        circuits.append(
            generate_tfim_circ_custom_graph(
                steps, h, J, dt, psi0, G, meas_basis, mirror
            )
        )
    return circuits

Estimar las magnetizaciones de un solo cúbit $\langle Z_i \rangle$

Para estudiar la dinámica del modelo, queremos medir la magnetización de cada cúbit, definida por el valor esperado $\langle Z_i \rangle = \langle \psi | Z_i | \psi \rangle$.

En simulaciones, podemos calcular esto directamente a partir de los resultados de medición. La función z_expectation procesa los conteos de cadenas de bits y devuelve el valor de $\langle Z_i \rangle$ para un índice de cúbit elegido. En hardware real, evaluamos la misma cantidad especificando el operador de Pauli utilizando la función generate_z_observables, y luego el backend calcula el valor esperado.

def z_expectation(counts, index):
    """
    counts: Dict of mitigated bitstrings.
    index: Index i in the single operator expectation value < II...Z_i...I > to be calculated.
    return:  < Z_i >
    """
    z_exp = 0
    tot = 0
    for bitstring, value in counts.items():
        bit = int(bitstring[index])
        sign = 1
        if bit % 2 == 1:
            sign = -1
        z_exp += sign * value
        tot += value

    return z_exp / tot

def generate_z_observables(nq):
    observables = []
    for i in range(nq):
        pauli_string = "".join(["Z" if j == i else "I" for j in range(nq)])
        observables.append(SparsePauliOp(pauli_string))
    return observables

observables = generate_z_observables(n_qubits)

Ahora definimos los parámetros para generar los circuitos trotterizados. En este tutorial, la red es una cadena de 20 triángulos conectados, lo que corresponde a un sistema de 41 cúbits.

all_circs_mirror = []
for num_triangles in [n_triangles]:
    for meas_basis in ["Z"]:
        A = connected_triangles_adj_matrix(num_triangles)
        G = nx.from_numpy_array(A)
        nq = len(G)
        d_ind_tot = 22
        dt = 2 * np.pi * 1 / 30 * 0.25
        J = 1
        h = -7
        all_circs_mirror.extend(
            trotter_circuits(G, d_ind_tot, J, h, dt, meas_basis, True)
        )
circs = all_circs_mirror

Paso 3: Ejecutar utilizando primitivas de Qiskit

Ejecutar simulación MPS

La lista de circuitos trotterizados se ejecuta utilizando el simulador matrix_product_state con una elección arbitraria de $4096$ disparos. El método MPS proporciona una aproximación eficiente de la dinámica del circuito, con una precisión determinada por la dimensión de enlaza elegida. Para los tamaños de sistema considerados aquí, la dimensión de enlaza predeterminada es suficiente para capturar la dinámica de magnetización con alta fidelidad. Los conteos en bruto se normalizan, y a partir de estos calculamos los valores esperados de un solo cúbit $\langle Z_i \rangle$ en cada paso de Trotter. Finalmente, calculamos el promedio sobre todos los cúbits para obtener una curva única que muestra cómo cambia la magnetización a lo largo del tiempo.

backend_sim = AerSimulator(method="matrix_product_state")

def normalize_counts(counts_list, shots):
    new_counts_list = []
    for counts in counts_list:
        a = {k: v / shots for k, v in counts.items()}
        new_counts_list.append(a)
    return new_counts_list

def run_sim(circ_list):
    shots = 4096
    res = backend_sim.run(circ_list, shots=shots)
    normed = normalize_counts(res.result().get_counts(), shots)
    return normed

sim_counts = run_sim(circs)

Ejecutar en hardware

service = QiskitRuntimeService()
backend = service.backend("ibm_marrakesh")

def run_qiskit(circ_list):
    shots = 4096
    pm = generate_preset_pass_manager(backend=backend)
    isa_circuits = [pm.run(qc) for qc in circ_list]
    sampler = Sampler(mode=backend)
    res = sampler.run(isa_circuits, shots=shots)
    res = [r.data.c.get_counts() for r in res.result()]
    normed = normalize_counts(res, shots)
    return normed

qiskit_counts = run_qiskit(circs)

Ejecutar en hardware con Fire Opal

Evaluamos la dinámica de magnetización en hardware cuántico real. Fire Opal proporciona una función de Qiskit que extiende la primitiva estándar Estimator de Qiskit Runtime con supresión automática de errores y gestión de rendimiento. Enviamos los circuitos trotterizados directamente a un backend de IBM® mientras Fire Opal se encarga de la ejecución con reconocimiento de ruido.

Preparamos una lista de pubs, donde cada elemento contiene un circuito y los observables de Pauli-Z correspondientes. Estos se pasan a la función estimator de Fire Opal, que devuelve los valores esperados $\langle Z_i \rangle$ para cada cúbit en cada paso de Trotter. Los resultados pueden luego promediarse sobre los cúbits para obtener la curva de magnetización del hardware.

backend_name = "ibm_marrakesh"
estimator_pubs = [(qc, observables) for qc in all_circs_mirror[:]]

# Run the circuit using the estimator
qctrl_estimator_job = perf_mgmt.run(
    primitive="estimator",
    pubs=estimator_pubs,
    backend_name=backend_name,
    options={"default_shots": 4096},
)

result_qctrl = qctrl_estimator_job.result()

Paso 4: Posprocesar y devolver el resultado en el formato clásico deseado

Finalmente, comparamos la curva de magnetización del simulador con los resultados obtenidos en hardware real. Graficar ambos lado a lado muestra cuán estrechamente la ejecución en hardware con Fire Opal coincide con la referencia sin ruido a lo largo de los pasos de Trotter.

def make_correlators(test_counts, nq, d_ind_tot):
    mz = np.empty((nq, d_ind_tot))
    for d_ind in range(d_ind_tot):
        counts = test_counts[d_ind]
        for i in range(nq):
            mz[i, d_ind] = z_expectation(counts, i)
    average_z = np.mean(mz, axis=0)
    return np.concatenate((np.array([1]), average_z), axis=0)

sim_exp = make_correlators(sim_counts[0:22], nq=nq, d_ind_tot=22)
qiskit_exp = make_correlators(qiskit_counts[0:22], nq=nq, d_ind_tot=22)

qctrl_exp = [ev.data.evs for ev in result_qctrl[:]]
qctrl_exp_mean = np.concatenate(
    (np.array([1]), np.mean(qctrl_exp, axis=1)), axis=0
)

def make_expectations_plot(
    sim_z,
    depths,
    exp_qctrl=None,
    exp_qctrl_error=None,
    exp_qiskit=None,
    exp_qiskit_error=None,
    plot_from=0,
    plot_upto=23,
):
    import numpy as np
    import matplotlib.pyplot as plt

    depth_ticks = [0, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22]

    d = np.asarray(depths)[plot_from:plot_upto]
    sim = np.asarray(sim_z)[plot_from:plot_upto]

    qk = (
        None
        if exp_qiskit is None
        else np.asarray(exp_qiskit)[plot_from:plot_upto]
    )
    qc = (
        None
        if exp_qctrl is None
        else np.asarray(exp_qctrl)[plot_from:plot_upto]
    )

    qk_err = (
        None
        if exp_qiskit_error is None
        else np.asarray(exp_qiskit_error)[plot_from:plot_upto]
    )
    qc_err = (
        None
        if exp_qctrl_error is None
        else np.asarray(exp_qctrl_error)[plot_from:plot_upto]
    )

    # ---- helper(s) ----
    def rmse(a, b):
        if a is None or b is None:
            return None
        a = np.asarray(a, dtype=float)
        b = np.asarray(b, dtype=float)
        mask = np.isfinite(a) & np.isfinite(b)
        if not np.any(mask):
            return None
        diff = a[mask] - b[mask]
        return float(np.sqrt(np.mean(diff**2)))

    def plot_panel(ax, method_y, method_err, color, label, band_color=None):
        # Noiseless reference
        ax.plot(d, sim, color="grey", label="Noiseless simulation")

        # Method line + band
        if method_y is not None:
            ax.plot(d, method_y, color=color, label=label)
            if method_err is not None:
                lo = np.clip(method_y - method_err, -1.05, 1.05)
                hi = np.clip(method_y + method_err, -1.05, 1.05)
                ax.fill_between(
                    d,
                    lo,
                    hi,
                    alpha=0.18,
                    color=band_color if band_color else color,
                    label=f"{label} ± error",
                )
        else:
            ax.text(
                0.5,
                0.5,
                "No data",
                transform=ax.transAxes,
                ha="center",
                va="center",
                fontsize=10,
                color="0.4",
            )

        # RMSE box (vs sim)
        r = rmse(method_y, sim)
        if r is not None:
            ax.text(
                0.98,
                0.02,
                f"RMSE: {r:.4f}",
                transform=ax.transAxes,
                va="bottom",
                ha="right",
                fontsize=8,
                bbox=dict(
                    boxstyle="round,pad=0.35", fc="white", ec="0.7", alpha=0.9
                ),
            )
        # Axes
        ax.set_xticks(depth_ticks)
        ax.set_ylim(-1.05, 1.05)
        ax.grid(True, which="both", linewidth=0.4, alpha=0.4)
        ax.set_axisbelow(True)
        ax.legend(prop={"size": 8}, loc="best")

    fig, axes = plt.subplots(1, 2, figsize=(10, 4), dpi=300, sharey=True)

    axes[0].set_title("Fire Opal (Q-CTRL)", fontsize=10)
    plot_panel(
        axes[0],
        qc,
        qc_err,
        color="#680CE9",
        label="Fire Opal",
        band_color="#680CE9",
    )
    axes[0].set_xlabel("Trotter step")
    axes[0].set_ylabel(r"$\langle Z \rangle$")
    axes[1].set_title("Qiskit", fontsize=10)
    plot_panel(
        axes[1], qk, qk_err, color="blue", label="Qiskit", band_color="blue"
    )
    axes[1].set_xlabel("Trotter step")

    plt.tight_layout()
    plt.show()

depths = list(range(d_ind_tot + 1))
errors = np.abs(np.array(qctrl_exp_mean) - np.array(sim_exp))

errors_qiskit = np.abs(np.array(qiskit_exp) - np.array(sim_exp))

make_expectations_plot(
    sim_exp,
    depths,
    exp_qctrl=qctrl_exp_mean,
    exp_qctrl_error=errors,
    exp_qiskit=qiskit_exp,
    exp_qiskit_error=errors_qiskit,
)

Referencias

[1] Graph coloring. Wikipedia. Retrieved September 15, 2025, from https://en.wikipedia.org/wiki/Graph_coloring

Encuesta del tutorial

Por favor, tómese un minuto para proporcionar comentarios sobre este tutorial. Sus opiniones nos ayudarán a mejorar nuestras ofertas de contenido y la experiencia del usuario.

Enlaza a la encuesta