SQD y SKQD

En este capítulo exploramos cómo los computadores cuánticos y clásicos trabajan juntos para abordar uno de los desafíos más importantes de la ciencia: la estimación precisa de la energía de moléculas y materiales.

Iskandar Sitdikov describe el enfoque algorítmico en el siguiente video.

Hamiltoniano

La clave de este problema es un operador matemático —el hamiltoniano, que describe la energía total de un sistema. Para fines computacionales, se puede pensar en este hamiltoniano como una gran matriz. Las soluciones buscadas —en particular, el estado fundamental del sistema— son los autovalores más bajos de esta matriz. Sin embargo, el desafío es que esta matriz hamiltoniana es muy grande en problemas prácticos. Crece exponencialmente con el tamaño del sistema y rápidamente se vuelve demasiado grande ($2^n$, donde $n$ es el número de qubits) para almacenarse o resolverse directamente incluso en los supercomputadores más potentes.

$$H = \begin{pmatrix} H_{0,0} & H_{0,1} & \cdots & H_{0,N-1} \\ H_{1,0} & H_{1,1} & \cdots & H_{1,N-1} \\ \vdots & \vdots & \ddots & \vdots \\ H_{N-1,0} & H_{N-1,1} & \cdots & H_{N-1,N-1} \end{pmatrix} \quad (N=2^n)$$

Para sortear este problema, utilizamos una estrategia poderosa conocida como el método del subespacio. En lugar de abordar la matriz completa, seleccionamos inteligentemente una porción pequeña y relevante —un "subespacio"— del cual asumimos que contiene la información más importante sobre la solución de baja energía buscada.

$$\underset{\text{Full Hamiltonian}}{\begin{pmatrix} \ddots & \vdots \\ \cdots & H_{i,j} & \cdots \\ & \vdots & \ddots \end{pmatrix}} \quad \xrightarrow{\text{Project}} \quad \underset{\text{Projected Hamiltonian}}{\tilde{H}} = \begin{pmatrix} \langle b_1 | H | b_1 \rangle & \cdots & \langle b_1 | H | b_L \rangle \\ \vdots & \ddots & \vdots \\ \langle b_L | H | b_1 \rangle & \cdots & \langle b_L | H | b_L \rangle \end{pmatrix} \quad \xrightarrow{\text{Diagonalize}} \quad \underset{\text{Eigenvalues}}{\begin{pmatrix} E_0 & & 0 \\ & \ddots & \\ 0 & & E_{L-1} \end{pmatrix}}$$

Una vez que este pequeño subespacio está definido por un conjunto de estados base $\{|b_i\rangle\}$, el hamiltoniano completo se proyecta sobre él para producir la nueva matriz más pequeña $\tilde{H}$. Cada elemento de esta matriz se calcula a partir de los estados base del subespacio y el hamiltoniano original como $\langle b_i | H | b_j \rangle$. Esta pequeña matriz puede entonces diagonalizarse fácilmente en un computador clásico, y los autovalores resultantes son nuestras energías estimadas.

Como puedes imaginar, el éxito de todo este enfoque depende fuertemente de la elección de un "buen" subespacio. Si nuestro subespacio no representa bien el verdadero estado fundamental, nuestro resultado final será incorrecto. Aquí es donde entran los computadores cuánticos: nos permiten preparar y muestrear estados cuánticos complejos diseñados para identificar estos subespacios importantes. Para problemas muy grandes, como estructuras químicas complejas o sitios de enlace, incluso la matriz proyectada puede ser difícil de diagonalizar. Tales problemas son por lo tanto ideales para combinar las fortalezas de los recursos de cómputo cuánticos y clásicos.

En las siguientes secciones, exploraremos dos algoritmos avanzados —SQD y SKQD— que utilizan la mecánica cuántica para encontrar y construir estos subespacios. Para una exploración más profunda, hay un curso completo en IBM Quantum Learning dedicado a estos temas en detalle. En el marco de este curso, mantenemos la explicación a un nivel general.

Diagonalización Cuántica Basada en Muestras

La Diagonalización Cuántica Basada en Muestras (SQD) es un potente algoritmo variacional que implementa el método del subespacio de manera cuántica. Evita procedimientos costosos y complejos como los tests de Hadamard utilizando un computador cuántico para preparar un estado de prueba y muestrear cadenas de bits que definen el subespacio para la diagonalización clásica.

El algoritmo SQD se puede dividir en los siguientes pasos:

Paso 1: Preparar el estado ansatz

Sea $H = \sum_{j=1}^Q \alpha_j P_j$ el hamiltoniano sobre $n$ qubits. Aunque el verdadero estado fundamental puede estar distribuido sobre todos los $2^n$ estados base, SQD es más efectivo en casos donde el estado fundamental puede ser bien aproximado por un subespacio disperso (un conjunto polinomialmente grande de cadenas de bits).

Para construir este subespacio, comenzamos con un estado de entrada $|\phi_0\rangle$, por ejemplo, el estado de Hartree-Fock (HF) en química. Luego aplicamos un circuito cuántico parametrizado $U(\theta)$, denominado ansatz.

Este diagrama ilustra el objetivo de un buen ansatz. El ansatz prepara un estado cuántico cuyo soporte (el conjunto de estados base de los que está compuesto) idealmente debería tener una gran superposición con el soporte del verdadero estado fundamental. Este circuito nos permite proyectar rápidamente el ansatz sobre vectores de la base computacional que se utilizarán en la diagonalización clásica. Dicho de otro modo: no necesitamos adivinar un ansatz que corresponda al estado fundamental; solo necesitamos que contenga los mismos estados base. La diagonalización clásica del hamiltoniano proyectado nos proporcionará entonces la superposición de estados base que mejor aproxima el estado fundamental.

Paso 2: Muestrear el subespacio

Al muestrear del circuito preparado por el ansatz, obtenemos una colección de cadenas de bits $\{b_j\}_{j=1}^L$. Estas cadenas de bits definen la base de nuestro subespacio elegido. El tiempo de ejecución cuántica para este paso está determinado por la profundidad del circuito y el número de muestras tomadas.

Paso 3: Proyectar y diagonalizar clásicamente

Usando las cadenas de bits muestreadas, proyectamos el hamiltoniano en el subespacio generado por ellas. Para cada par de cadenas de bits $(j, k)$, calculamos clásicamente el elemento de matriz $\tilde{H}_{jk} = \langle b_j | H | b_k \rangle$. Dado que los operadores de Pauli son dispersos, este paso es clásicamente eficiente para hamiltonianos físicos. La pequeña matriz resultante $\tilde{H}$ se diagonaliza entonces en un procesador clásico para estimar el estado fundamental y su energía.

Paso 4: Optimizar el ansatz (opcional)

El proceso puede hacerse iterativo. Tratando la energía estimada del estado fundamental como función de costo, podemos optimizar los parámetros del circuito ($\theta$) mediante métodos como el descenso de gradiente para mejorar el ansatz y, por lo tanto, la aproximación de la energía en la siguiente iteración.

Ventajas clave de SQD

SQD ofrece varias propiedades poderosas que lo convierten en un candidato líder para la demostración de la ventaja cuántica:

• Fuerte robustez frente al ruido: Supongamos que el verdadero estado fundamental se apoya solo en dos cadenas de bits. Si estas se muestrean de alguna manera —incluso si su superposición con nuestro ansatz es pequeña—, la diagonalización les asigna las ponderaciones correspondientes e ignora efectivamente todas las demás cadenas de bits superfluas y ruidosas que también podrían haberse muestreado. Este filtrado inherente hace que SQD sea particularmente robusto frente al ruido.
• Verificabilidad clásica: A diferencia de QPE o VQAs, SQD proporciona una aproximación clásica al estado fundamental. Esto significa que cualquier persona que conozca la lista de cadenas de bits y sus ponderaciones puede recalcular y verificar directamente la estimación de energía en un computador clásico.

SQD ya se ha utilizado para estimar la energía de disociación del estado fundamental del N$_2$ y las propiedades electrónicas de clústeres [2Fe-2S] y [4Fe-4S] [2], con circuitos de hasta 77 qubits y 10.570 puertas.

Preguntas de comprensión

Verdadero o Falso: SQD se puede aplicar a sistemas químicos.

Respuesta:

Verdadero

Preguntas de comprensión

Llama al conjunto de todos los vectores de la base computacional de tu ansatz $A$. Llama al conjunto de todos los vectores de la base computacional del verdadero estado fundamental de tu sistema $G$. ¿Cuál de las siguientes afirmaciones corresponde a un "buen" ansatz? Selecciona todas las que apliquen.

(a) $A \subset G \\$ (b) $A \subseteq G\\$ (c) $G \subset A\\$ (d) $G \subseteq A\\$

Respuesta:

(c) y (d)

SKQD (Diagonalización Cuántica de Krylov Basada en Muestras)

La Diagonalización Cuántica de Krylov Basada en Muestras (SKQD) es otro potente algoritmo cuántico basado en muestreo que se basa en los principios de SQD. Aunque su objetivo es el mismo —encontrar un buen subespacio para la diagonalización—, SKQD utiliza un método más estructurado para generar las cadenas de bits, particularmente para problemas como hamiltonianos de red.

La idea central de SKQD es que, en lugar de optimizar un circuito parametrizado para un buen ansatz, se puede lograr una convergencia demostrable al estado fundamental muestreando de un conjunto de estados generados por la evolución temporal natural del sistema mismo: el subespacio de Krylov. El algoritmo SKQD se puede dividir en los siguientes pasos:

Paso 1: Construir el subespacio de Krylov con evolución temporal

El proceso comienza con un estado inicial $|\phi_0\rangle$. Es importante señalar que este estado inicial no necesita tener una "buena" superposición con el estado fundamental. Solo necesita ser "polinomialmente grande", es decir, describible por un polinomio en el tamaño del sistema. El algoritmo mismo impulsa el estado cada vez más cerca del estado fundamental del sistema. SKQD aplica el operador de evolución temporal $e^{-iHt}$ para diferentes duraciones de tiempo. Esto genera un conjunto de $d$ estados cuánticos diferentes, definidos como:

$$|\phi_j\rangle = e^{-i \,\delta t j H}|\phi_0\rangle, \quad \text{para } j = 0, 1, \dots, d-1 \quad \text{}$$

Esta colección de estados evolucionados temporalmente forma una base de Krylov. Este paso es particularmente efectivo para hamiltonianos de red, donde el número de términos en el hamiltoniano no es grande. Para problemas de química, esta evolución temporal puede llevar a circuitos muy profundos, por lo que para tales casos a menudo se recomienda SQD.

Paso 2: Muestrear de los estados base de Krylov

A continuación, se recopilan muestras de cadenas de bits de cada uno de los $d$ estados diferentes ($|\phi_0\rangle, |\phi_1\rangle, \dots, |\phi_{d-1}\rangle$) preparados en el paso anterior. Todas estas cadenas de bits se combinan para formar la base del subespacio.

Paso 3: Proyectar y diagonalizar clásicamente

Este paso es idéntico al paso correspondiente en SQD. Las cadenas de bits recopiladas se utilizan para proyectar el hamiltoniano completo en el subespacio generado por ellas. La pequeña matriz resultante $\tilde{H}$ se diagonaliza entonces en un computador clásico para encontrar la energía del estado fundamental.

Ventajas y garantías clave de SKQD

El enfoque estructurado de SKQD ofrece ventajas únicas:

• Convergencia demostrable: La principal ventaja de SKQD es su garantía teórica de convergencia bajo condiciones específicas y claramente definidas. Si el verdadero estado fundamental es disperso (puede ser bien aproximado por un número polinomial de cadenas de bits) y la brecha energética al primer estado excitado no es demasiado pequeña, el método es demostrablemente eficiente. Bajo estas condiciones, SKQD garantiza que encontrará las cadenas de bits cruciales de las que se compone el estado fundamental, y puede aproximar la energía del estado fundamental con alta precisión. Para ello, solo se requiere un número polinomial de experimentos cuánticos y disparos. Esta garantía pone el enfoque basado en muestreo sobre una base teórica rigurosa, similar a métodos establecidos como la estimación de fase cuántica.

• Ventajas compartidas con SQD: Al igual que SQD, SKQD también posee la propiedad de robustez frente al ruido. Es decir: siempre que todas las buenas cadenas de bits estén presentes en el conjunto de cadenas de bits muestreadas, la diagonalización asigna peso casi nulo a las cadenas de bits incorrectas, lo que hace al procedimiento robusto frente al ruido. Dado que la solución es producida por un HPC clásico, la energía de la solución es verificable clásicamente.

En experimentos, SKQD se ha utilizado con hasta 70 qubits y miles de puertas para investigar el estado fundamental de modelos complejos de Anderson de 4 impurezas, logrando una excelente concordancia con métodos clásicos de vanguardia como DMRG.[1]

Preguntas de comprensión

¿Qué parte del algoritmo SKQD lo hace más adecuado para problemas físicos como redes de espín que para problemas químicos? ¿Por qué?

Respuesta:

La evolución temporal requiere circuitos de Trotter, que son muy profundos para hamiltonianos complicados y no dispersos. Las interacciones de redes de espín se describen mediante matrices de espín que corresponden a las matrices de Pauli. Por lo tanto, los hamiltonianos para redes de espín se expresan de forma más compacta en matrices de Pauli, especialmente para interacciones de vecinos más cercanos.

SQD y SKQD como computación heterogénea

Para unir todo, podemos representar los algoritmos basados en muestreo como una combinación de diferentes modelos de programación sobre un conjunto de recursos heterogéneos. Podemos representar nuestro algoritmo, por ejemplo, como un flujo de trabajo de tareas.

Esta figura ilustra el flujo de trabajo fundamental de cuatro etapas. Primero, hay una tarea de preparación del circuito cuántico que se superpone con nuestro estado objetivo, seguida de una tarea de transpilación que solo requiere recursos clásicos. A continuación, hay una tarea que utiliza primitivos para ejecutar nuestro circuito cuántico, que requiere recursos cuánticos. Finalmente, hay una tarea de posprocesamiento que puede ser en sí misma un algoritmo de diagonalización paralelo ejecutándose en múltiples nodos.

Además, puede ser útil ejecutar uno de estos algoritmos múltiples veces mientras variamos nuestro ansatz, o ejecutarlos completamente en paralelo con diferentes poblaciones.

Como se muestra arriba, puedes ejecutar múltiples flujos de trabajo simultáneamente haciendo lo siguiente:

• Variar los parámetros o la estructura del ansatz para encontrar el más efectivo.
• Comenzar con diferentes estados iniciales o configuraciones ("poblaciones") para evitar mínimos locales y asegurar un resultado más robusto.

Este enfoque multicapa, que combina heterogeneidad basada en tareas con paralelismo a nivel de flujo de trabajo, es la clave para aprovechar al máximo el potencial de estos algoritmos.

Práctica de programación

Practiquemos ahora el algoritmo SKQD y demostremos el flujo de trabajo heterogéneo descrito anteriormente. El proceso se divide en cuatro etapas distintas, cada una con su propio script Python y un script de shell asociado para el envío de trabajos.

Mapping (mapping.py y mapping.sh)

El primer paso en nuestro flujo de trabajo consiste en definir el problema físico y mapearlo a un conjunto de circuitos cuánticos.

mapping.py define los parámetros para un problema físico específico —en este caso, un modelo de impurezas de Anderson con siete sitios de baño (n_bath = 7). Construye las integrales de un cuerpo (h1e) y de dos cuerpos (h2e) que describen el hamiltoniano del sistema.

...

n_bath = 7  # number of bath sites

...

# One body matrix elements in the "position" basis
h1e = -t * np.diag(np.ones(n_bath), k=1) - t * np.diag(np.ones(n_bath), k=-1)
h1e[impurity_index, impurity_index + 1] = -V
h1e[impurity_index + 1, impurity_index] = -V
h1e[impurity_index, impurity_index] = eps

# Two body matrix elements in the "position" basis
h2e = np.zeros((n_bath + 1, n_bath + 1, n_bath + 1, n_bath + 1))
h2e[impurity_index, impurity_index, impurity_index, impurity_index] = U

...

# The one-body time evolution
free_fermion_evolution = ffsim.qiskit.OrbitalRotationJW(n_modes, Utar)

# The two-body time evolution
def append_diagonal_evolution(dt, U, impurity_qubit, num_orb, q_circuit):
    """Append two-body time evolution to a quantum circuit."""
    if U != 0:
        q_circuit.append(
            CPhaseGate(-dt / 2 * U),
            [impurity_qubit, impurity_qubit + num_orb],
        )

A continuación, el script genera los circuitos cuánticos necesarios para el algoritmo SKQD. Comienza creando un estado inicial (initial_state) y luego aplica operadores de evolución temporal para un número variable de pasos (d = 8) para producir los diferentes estados base de Krylov $|\phi_j\rangle = (e^{-iHt})^j |\phi_0\rangle$.

# The reference state
def initial_state(q_circuit, norb, nocc):
    """Prepare an initial state."""
    for i in range(nocc):
        q_circuit.append(XGate(), [i])
        q_circuit.append(XGate(), [norb + i])
    rot = XXPlusYYGate(np.pi / 2, -np.pi / 2)

    for i in range(3):
        for j in range(nocc - i - 1, nocc + i, 2):
            q_circuit.append(rot, [j, j + 1])
            q_circuit.append(rot, [norb + j, norb + j + 1])
    q_circuit.append(rot, [j + 1, j + 2])
    q_circuit.append(rot, [norb + j + 1, norb + j + 2])

...

# Generate the initial state
qubits = QuantumRegister(2 * n_modes, name="q")
init_state = QuantumCircuit(qubits)
initial_state(init_state, n_modes, n_modes // 2)

...

d = 8  # Number of Krylov basis states
circuits = []
for i in range(d):
    circ = init_state.copy()
    circuits.append(circ)
    for _ in range(i):
        append_diagonal_evolution(dt, U, impurity_index, n_modes, circ)
        circ.append(free_fermion_evolution, qubits)
        append_diagonal_evolution(dt, U, impurity_index, n_modes, circ)
    circ.measure_all()

print(circuits[0].draw(scale=0.4, fold=-1))

El script guarda la lista de 8 circuitos generados (cada uno con mediciones adjuntas) en un archivo llamado circuits.qpy.

mapping.sh es un script por lotes de Slurm utilizado para enviar el trabajo mapping.py. Dado que se trata de un cálculo clásico, se solicitan recursos de una partición CPU estándar (--partition=normal).

#!/bin/bash
#
#SBATCH --job-name=sqd-mapping
#SBATCH --output=sqd-mapping.out
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=1
#SBATCH --cpus-per-task=1
#SBATCH --partition=normal

srun python /data/ch4/sqd/mapping.py

Optimización (optimization.py y optimization.sh)

Una vez que tenemos nuestros circuitos, deben optimizarse y compilarse para que puedan ejecutarse eficientemente en el hardware cuántico objetivo.

En optimization.py, el script primero carga el archivo circuits.qpy creado en la etapa de mapping y obtiene información de recursos cuánticos a través de QRMI(), un gestor de recursos cuánticos. Luego utiliza generate_preset_pass_manager de Qiskit con un alto nivel de optimización (optimization_level=3) para convertir los circuitos abstractos y lógicos en circuitos ISA (Instruction Set Architecture). Este proceso reescribe los circuitos utilizando las puertas nativas del hardware y los optimiza para reducir la profundidad y minimizar los errores.

...
qrmi = QRMI()
resources = qrmi.resources()
quantum_resource = resources[0]
target = quantum_resource.target

pass_manager = generate_preset_pass_manager(
    optimization_level=3,
    target=target
)
isa_circuits = pass_manager.run(circuits)

Los circuitos transpilados y listos para el hardware se guardan en un nuevo archivo isa_circuits.qpy.

Al igual que el script de mapping, este trabajo de Slurm también se ejecuta en una partición CPU clásica (--partition=normal), ya que la transpilación es una tarea clásica.

#!/bin/bash
#
#SBATCH --job-name=sqd-mapping
#SBATCH --output=sqd-mapping.out
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=1
#SBATCH --cpus-per-task=1
#SBATCH --partition=normal

srun python /data/ch4/sqd/mapping.py

Ejecución (execution.py y execution.sh)

Esta es la única etapa en la que se utiliza un computador cuántico. Aquí se ejecutan los circuitos optimizados y se recopilan muestras de medición.

execution.py carga el archivo optimizado isa_circuits.qpy e inicializa un primitivo SamplerV2 conectado a un recurso cuántico. Luego llama a sampler.run() para ejecutar los circuitos en la QPU con un número especificado de disparos (shots=500).

...

qrmi = QRMI()
resources = qrmi.resources()
quantum_resource = resources[0]

# Sample from the circuits
noisy_sampler = Sampler(quantum_resource)
job = noisy_sampler.run(isa_circuits, shots=500)

Al final de la ejecución, los resultados medidos (cadenas de bits) de todos los circuitos se recopilan y combinan, y sus conteos se guardan en un archivo counts.json.

El script de Slurm execution.sh difiere de los demás en esta etapa. Solicita ejecutarse en la partición cuántica (--partition=quantum) y solicita específicamente una QPU (--gres=qpu:1).

#!/bin/bash
#
#SBATCH --job-name=sqd-execution
#SBATCH --output=sqd-execution.out
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=1
#SBATCH --cpus-per-task=1
#SBATCH --partition=quantum
#SBATCH --gres=qpu:1

srun python /data/ch4/sqd/execution.py

Posprocesamiento (postprocessing.py y postprocessing.sh)

En el último paso, volvemos al computador clásico para analizar los datos del experimento cuántico y calcular el resultado final: la energía del estado fundamental de nuestro sistema objetivo.

postprocessing.py primero lee el archivo counts.json con los resultados de las mediciones. Luego reconstruye el hamiltoniano del modelo de Anderson (con los mismos parámetros que en mapping.py). A continuación, pasa las cadenas de bits medidas y la definición del hamiltoniano a la función diagonalize_fermionic_hamiltonian. Esta función ejecuta la lógica central de SKQD: utiliza las cadenas de bits para construir el hamiltoniano proyectado $\tilde{H}$ y lo diagonaliza para encontrar la energía del estado fundamental.

...

def callback(results: list[SCIResult]):
    result_history.append(results)
    iteration = len(result_history)
    print(f"Iteration {iteration}")
    for i, result in enumerate(results):
        print(f"\tSubsample {i}")
        print(f"\t\tEnergy: {result.energy}")
        print(f"\t\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}")

rng = np.random.default_rng(24)
result = diagonalize_fermionic_hamiltonian(
    h1e,
    h2e,
    bit_array,
    samples_per_batch=300,
    norb=n_modes,
    nelec=nelec,
    num_batches=3,
    max_iterations=10,
    symmetrize_spin=True,
    callback=callback,
    seed=rng,
)

Finalmente, muestra la energía SKQD calculada y la compara con la energía exacta conocida para este problema, mostrando así el error absoluto final del cálculo.

El script de trabajo final se ejecuta en una partición clásica (--partition=normal), ya que todo el análisis es de naturaleza clásica. Para subespacios grandes, este paso puede requerir más recursos HPC clásicos.

#!/bin/bash
#
#SBATCH --job-name=sqd-postprocessing
#SBATCH --output=sqd-postprocessing.out
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=1
#SBATCH --cpus-per-task=1
#SBATCH --partition=normal

srun python /data/ch4/sqd/postprocessing.py

Resumen

¡Eso es todo! Ahora hemos trabajado a través de varios conceptos y ejemplos que pueden ayudarte a comenzar con la gestión de programas híbridos complejos. Por supuesto, esto es solo el comienzo de lo que puedes lograr con la combinación de recursos cuánticos y HPC clásicos.

Para explorar más casos de uso y algoritmos, visita nuestra documentación y tutoriales en IBM Quantum Platform, y asegúrate de revisar los recursos compartidos en la próxima lección para obtener más información sobre algoritmos y software para científicos computacionales y administradores de centros de datos.

Referencias

[1] Quantum-Centric Algorithm for Sample-Based Krylov Diagonalization. https://arxiv.org/abs/2501.09702

[2] Chemistry beyond the scale of exact diagonalization on a quantum-centric supercomputer. https://www.science.org/doi/10.1126/sciadv.adu9991