Resolva o modelo de Sherrington-Kirkpatrick com o ParityQC Parity Twine Optimizer

Estimativa de uso: 10 segundos em um processador Nighthawk r2. (OBSERVAÇÃO: Esta é apenas uma estimativa. Seu tempo de execução pode variar.)

Resultados de aprendizado

• Use o Parity Twine Optimizer para resolver o modelo de Sherrington-Kirkpatrick.

• Aprenda quais opções estão disponíveis para o Parity Twine Optimizer e quais resultados são gerados.

Contexto

Este tutorial demonstra como resolver o modelo de Sherrington-Kirkpatrick usando o ParityQC Parity Twine Optimizer.

Ele fornece código para formular o problema inteiramente de forma local em um formato que pode se conectar com o Parity Twine Optimizer.

O modelo de Sherrington-Kirkpatrick

O modelo de Sherrington-Kirkpatrick (SK) é um modelo fundamental na mecânica estatística, especificamente no estudo de vidros de spin. Ao contrário do modelo de Ising padrão, no qual as interações geralmente se limitam aos vizinhos mais próximos, o modelo SK é um modelo de alcance infinito, o que significa que cada spin interage com todos os outros spins do sistema. Isso leva a uma paisagem de energia altamente complexa e "acidentada", caracterizada por muitos mínimos locais, que é a marca registrada do comportamento vítreo.

A característica central do modelo SK reside na frustração. No modelo, as forças de interação $J_{ij}$ entre spins são distribuídas aleatoriamente entre valores positivos e negativos. Isso leva a situações (por exemplo, em arranjos triangulares) em que os spins não podem ser organizados para minimizar todas as interações simultaneamente. No modelo SK, como cada spin interage com todos os outros spins, essa frustração se agrava globalmente, levando a uma rede de restrições conflitantes.

Formulação matemática

O estado do sistema é definido por um conjunto de $N$ spins de Ising, $s_i \in \{+1, -1 \}$. A energia de uma configuração específica é dada pelo Hamiltoniano:

$H = - \sum_{1 \leq i \leq j \leq N} J_{ij} s_i s_j$

onde $J_{ij}$ é a força de acoplamento entre o spin $i$ e o spin $j$.

No modelo SK, os acoplamentos $J_{ij}$ são variáveis aleatórias independentes e identicamente distribuídas. Para garantir que a energia permaneça extensiva (proporcional a $N$) à medida que $N \rightarrow \infty$, a variância dos acoplamentos deve escalar com o número de partículas:

$J_{ij} \sim \mathcal{N} \left( 0, \frac{J^2}{N} \right)$

O estado fundamental para um determinado $H$ é a configuração específica de spins ($s_1, s_2, ...s_n$) que minimiza a energia. Encontrar o estado fundamental é um problema de otimização NP-difícil.

Requisitos

Antes de iniciar este tutorial, certifique-se de que o seguinte esteja instalado:

• Qiskit Functions Catalog IBM Client (pip install qiskit-ibm-catalog)

• Qiskit addon Optimization Mapper (pip install qiskit_addon_opt_mapper)

• NumPy (pip install numpy)

Você também precisa de permissão para acessar a função ParityQC Twine Optimizer. Para solicitar acesso, preencha este formulário.

Configuração

(Este código pressupõe que você já salvou sua conta no seu ambiente local.)

Primeiro, importe todos os pacotes necessários para este tutorial.

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

_needed = {"numpy": "numpy", "qiskit_ibm_catalog": "qiskit-ibm-catalog"}
_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")

import numpy as np

from qiskit_ibm_catalog import QiskitFunctionsCatalog

Carregue o Parity Twine Optimizer do catálogo do Qiskit Functions:

catalog = QiskitFunctionsCatalog(channel="ibm_quantum_platform")
function = catalog.load("parityqc/parity-twine-optimizer")

Etapa 1: Defina o problema como uma função objetivo

Em vez de obter o problema SK a partir de uma biblioteca, como fazemos para o problema Market Split, você o formula diretamente.

A função generate_sk_problem formula o problema SK diretamente no formato de dicionário necessário. A única entrada necessária é n, o número de spins no modelo.

def generate_sk_problem(
    n: int,
    coupling_mean: float = 0.0,
    coupling_std: float = 1.0,
    local_fields_mean: float = 0.0,
    local_fields_std: float = 0.0,
    edge_density: float = 1.0,
    ensure_extensivity: bool = False,
    seed: int | None = None,
) -> dict:
    """Generate the Sherrington-Kirkpatrick (SK) model with varying
     edge density.

    Samples couplings and local fields via ``generate_couplings_sk_model``
    and assembles the corresponding Ising Hamiltonian

        H = -∑_{i<j} J_ij z_i z_j - ∑_i h_i z_i,

    where z_i ∈ {-1, +1}.

    Args:
        n: Number of spins (>= 2).
        coupling_mean: Mean coupling before optional SK scaling.
        coupling_std: Coupling std before optional SK scaling.
        local_fields_mean: Mean longitudinal field.
        local_fields_std: Std of the longitudinal fields.
        edge_density: Fraction of non-zero couplings, in ``[2/n, 1]``.
        ensure_extensivity: Whether to apply the SK 1/n scaling.
        seed: random number generator seed.

    Returns:
        An encoding the SK Hamiltonian.

    Raises:
        ValueError: If ``n < 2``, ``coupling_std < 0``, ``local_fields_std < 0``,
            or ``edge_density`` is outside ``[2/n, 1]``.
    """
    couplings, local_fields = _generate_couplings_sk_model(
        n=n,
        coupling_mean=coupling_mean,
        coupling_std=coupling_std,
        local_fields_mean=local_fields_mean,
        local_fields_std=local_fields_std,
        edge_density=edge_density,
        ensure_extensivity=ensure_extensivity,
        seed=seed,
    )

    # Handle quadratic terms: coupling[i, j] * zj[i] * zj[j]
    # Only iterate over the upper triangle (i < j)
    sk_problem = {
        str((i, j)): float(couplings[i, j])
        for i in range(n)
        for j in range(i + 1, n)
        if couplings[i, j] != 0
    }

    # Handle linear terms: local_fields[i] * zj[i]
    sk_problem.update(
        {
            str((i,)): float(local_fields[i])
            for i in range(n)
            if local_fields[i] != 0
        }
    )

    return sk_problem

Use a função _generate_couplings_sk_model para calcular os termos de acoplamento aleatórios no modelo SK. Para ter mais controle sobre os acoplamentos, você pode usar argumentos opcionais, que são explicados na docstring da função.

def _generate_couplings_sk_model(
    n: int,
    coupling_mean: float = 0.0,
    coupling_std: float = 1.0,
    local_fields_mean: float = 0.0,
    local_fields_std: float = 0.0,
    edge_density: float = 1.0,
    ensure_extensivity: bool = False,
    seed: int | None = None,
) -> tuple[np.ndarray, np.ndarray]:
    """Generate random couplings and local fields for an Ising / SK model.

    Couplings are Gaussian. With ``ensure_extensivity=True`` they follow the
    Sherrington-Kirkpatrick scaling ``J_ij ~ N(coupling_mean/n, coupling_std^2/n)``
    (extensive energy, O(n)); otherwise ``J_ij ~ N(coupling_mean, coupling_std^2)``
    (energy O(n^2)). Fields are ``h_i ~ N(local_fields_mean, local_fields_std^2)``.

    ``edge_density`` sets the fraction of the ``n*(n-1)/2`` possible couplings that
    are non-zero (1 = fully dense). The kept edges always include a random spanning
    tree, so the interaction graph is guaranteed connected. This requires at least
    ``n-1`` edges, so ``edge_density`` must be at least ``2/n``.

    Args:
        n: Number of spins (>= 2).
        coupling_mean: Mean coupling before optional SK scaling.
        coupling_std: Coupling std before optional SK scaling.
        local_fields_mean: Mean longitudinal field.
        local_fields_std: Std of the longitudinal fields.
        edge_density: Fraction of non-zero couplings, in ``[2/n, 1]``.
        ensure_extensivity: Whether to apply the SK 1/n scaling.
        seed: random number generator seed.

    Returns:
        Tuple ``(couplings, fields)``: a symmetric ``(n, n)`` matrix with zero
        diagonal, and an ``(n,)`` field vector.

    Raises:
        ValueError: If ``n < 2``, ``coupling_std < 0``, ``local_fields_std < 0``,
            or ``edge_density`` is outside ``[2/n, 1]``.
    """
    if n < 2:
        raise ValueError(f"n must be >= 2, got {n}")
    if coupling_std < 0 or local_fields_std < 0:
        raise ValueError(
            "coupling_std and local_fields_std must be non-negative"
        )

    # A connected graph on n nodes needs at least n-1 of the n*(n-1)/2 possible
    # edges, so edge_density has a hard lower bound of 2/n.
    min_edge_density = 2.0 / n
    if not min_edge_density <= edge_density <= 1.0:
        raise ValueError(
            f"edge_density must be in [{min_edge_density:.4g}, 1] for n={n} "
            f"(at least n-1 edges are needed to keep the graph connected), "
            f"got {edge_density}"
        )

    rng = np.random.default_rng(seed)

    j_loc, j_scale = (
        (coupling_mean / n, coupling_std / np.sqrt(n))
        if ensure_extensivity
        else (coupling_mean, coupling_std)
    )

    upper_idx = np.triu_indices(n, k=1)
    n_edges = len(upper_idx[0])

    # Select which edges are present.
    if edge_density < 1.0:
        n_keep = int(round(edge_density * n_edges))
        # Map each (i, j) node pair to its position in the flat upper-triangle list.
        pair_to_flat = {
            (int(i), int(j)): idx
            for idx, (i, j) in enumerate(
                zip(upper_idx[0], upper_idx[1], strict=False)
            )
        }

        # Random spanning tree: node perm[k] links to a random earlier node.
        perm = rng.permutation(n)
        keep = np.zeros(n_edges, dtype=bool)
        for k in range(1, n):
            child, parent = perm[k], perm[rng.integers(0, k)]
            i, j = min(child, parent), max(child, parent)
            keep[pair_to_flat[(int(i), int(j))]] = True

        # Fill the remaining budget with random non-tree edges.
        remaining = n_keep - (n - 1)
        if remaining > 0:
            keep[
                rng.choice(
                    np.flatnonzero(~keep), size=remaining, replace=False
                )
            ] = True
    else:
        keep = np.ones(n_edges, dtype=bool)

    n_present = int(keep.sum())
    if j_scale == 0.0:
        vals = np.full(n_present, j_loc)
    else:
        vals = rng.normal(loc=j_loc, scale=j_scale, size=n_present)

    couplings = np.zeros((n, n))
    couplings[upper_idx[0][keep], upper_idx[1][keep]] = vals
    couplings += couplings.T  # symmetrize; diagonal stays zero

    fields = (
        np.full(n, local_fields_mean)
        if local_fields_std == 0.0
        else rng.normal(loc=local_fields_mean, scale=local_fields_std, size=n)
    )

    return couplings, fields

Etapa 2: Resolva o problema usando o Parity Twine Optimizer

Com as funções acima, você pode configurar o problema SK e encontrar uma solução usando o Twine Optimizer e um backend IBM Quantum® escolhido.

Para executar a função, escolha um backend adequado; por exemplo, ibm_phoenix.

Opcionalmente, use options para obter controle adicional sobre a submissão:

options = {
    "shots": 100000,
    "postprocessing_level": 1,
    "transpile_only": False,
    "job_tags": ["sk"],
}

onde shots é um número inteiro que especifica o número de execuções do circuito, postprocessing_level determina se o pós-processamento é aplicado ao resultado, transpile_only escolhe se o problema é apenas transpilado para um circuito (e não resolvido), e job_tags é a etiqueta para identificar o job na IBM Quantum Platform.

O tamanho do modelo SK é definido por $N$, o número de spins que interagem. Depois de escolher $N$, o código acima gera o problema para n_spins.

Execute o otimizador:

n_spins = 50
sk_problem = generate_sk_problem(n_spins)

function_job = function.run(
    problem=sk_problem,
    variable_type="spin",
    backend_name="ibm_phoenix",
    options=options,
)
print(f"Job ID: {function_job.job_id}")

Verifique o status do job:

# Monitor the job status
function_job.status()

Recupere os resultados:

result = function_job.result()

result

O resultado tem o seguinte formato:

{
    'solution': {'0': 1, '1': 1, '10': 1, '11': 1, ... },
    'objective_value':  -240.5425312543882,
    'solution_bitstring': '00001101110100100111001110101001101111111011001110',
    'metadata': {
        'circuit_metrics': {
            'depth': 523,
            'gate_count': 10118,
            'two_qubit_gate_depth': 196,
            'two_qubit_gate_count': 2499,
            'num_qubits': 50,
            'operations': {'sx': 3353, 'rz': 3320, 'cz': 2499, 'delay': 894, 'measure': 50, 'x': 2},
        },
        'solver_info': {
            'variable_mapping': {'0': 0, '1': 1, '10': 2, '11': 3, ... },
            'bitstring_distributions': {
                'before_postprocessing': {'011101110010110111001110011000': 1, ... },
                'after_postprocessing': {'011011110000110101001111011000': 1, ... }
            },
            'best_parameters': {
                'beta': [-0.4602084830507902],
                'gamma': [1.8500357096574955]
            }
        },
        'resource_usage': {
            'RUNNING: MAPPING': {'CPU_TIME': 290.272},
            'RUNNING: OPTIMIZING_FOR_HARDWARE': {'CPU_TIME': 0.494},
            'RUNNING: WAITING_FOR_QPU': {'CPU_TIME': 8.775},
            'RUNNING: EXECUTING_QPU': {'QPU_TIME': 31.0},
            'RUNNING: POST_PROCESSING': {'CPU_TIME': 162.96},
        },
    }
}

onde o dicionário solution corresponde aos qubits definidos no problema e fornece seus valores de spin otimizados para o Hamiltoniano do modelo SK. Essa sequência específica de spins na solução ótima representa a configuração que minimiza a energia total do sistema com base nas forças de interação aleatórias fornecidas. No modelo SK, isso pode ser considerado o estado de energia mais baixo de um sistema magnético desordenado.

metadata fornece informações sobre a transpilação (contagem/profundidade de portas de dois qubits, portas usadas, qubits ativos) e vários tempos de execução.

Próximos passos

Recomendações

• Solicite acesso à função preenchendo este formulário.
• Visite a referência da API para esta Qiskit Function.
• Leia o guia.
• Experimente o tutorial para aplicar o Parity Twine Optimizer ao problema Market Split.
• Revise o preprint do ArXiv Connectivity-aware Synthesis of Quantum Algorithms, Drier et al. (2025).