Löse das Sherrington-Kirkpatrick-Modell mit dem ParityQC Parity Twine Optimizer

Geschätzte Nutzungsdauer: 10 Sekunden auf einem Nighthawk r2-Prozessor. (HINWEIS: Dies ist nur eine Schätzung. Deine Laufzeit kann abweichen.)

Lernziele

• Verwende den Parity Twine Optimizer, um das Sherrington-Kirkpatrick-Modell zu lösen.

• Lerne, welche Optionen für den Parity Twine Optimizer verfügbar sind und welche Ergebnisse ausgegeben werden.

Hintergrund

Dieses Tutorial zeigt, wie du das Sherrington-Kirkpatrick-Modell mit dem ParityQC Parity Twine Optimizer löst.

Es stellt Code bereit, um das Problem vollständig lokal in einer Form zu formulieren, die mit dem Parity Twine Optimizer zusammenarbeiten kann.

Das Sherrington-Kirkpatrick-Modell

Das Sherrington-Kirkpatrick-Modell (SK-Modell) ist ein grundlegendes Modell der statistischen Mechanik, insbesondere im Bereich der Spinglas-Forschung. Anders als beim Standard-Ising-Modell, bei dem Wechselwirkungen üblicherweise auf die nächsten Nachbarn beschränkt sind, handelt es sich beim SK-Modell um ein Modell mit unendlicher Reichweite, das heißt, jeder Spin wechselwirkt mit jedem anderen Spin im System. Dies führt zu einer hochkomplexen, „zerklüfteten“ Energielandschaft, die durch viele lokale Minima gekennzeichnet ist – das Kennzeichen von glasartigem Verhalten.

Das Kernmerkmal des SK-Modells liegt in der Frustration. Im Modell sind die Wechselwirkungsstärken $J_{ij}$ zwischen Spins zufällig zwischen positiven und negativen Werten verteilt. Dies führt zu Situationen (zum Beispiel bei dreieckigen Anordnungen), in denen Spins nicht so angeordnet werden können, dass alle Wechselwirkungen gleichzeitig minimiert werden. Da im SK-Modell jeder Spin mit jedem anderen Spin wechselwirkt, wird diese Frustration global verstärkt, was zu einem Geflecht widersprüchlicher Bedingungen führt.

Mathematische Formulierung

Der Zustand des Systems wird durch eine Menge von $N$ Ising-Spins definiert, $s_i \in \{+1, -1 \}$. Die Energie einer bestimmten Konfiguration ist durch den Hamiltonoperator gegeben:

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

wobei $J_{ij}$ die Kopplungsstärke zwischen Spin $i$ und Spin $j$ ist.

Im SK-Modell sind die Kopplungen $J_{ij}$ unabhängige und identisch verteilte Zufallsvariablen. Um sicherzustellen, dass die Energie extensiv (proportional zu $N$) bleibt, wenn $N \rightarrow \infty$, muss die Varianz der Kopplungen mit der Anzahl der Teilchen skalieren:

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

Der Grundzustand für ein gegebenes $H$ ist die spezifische Konfiguration von Spins ($s_1, s_2, ...s_n$), die die Energie minimiert. Das Finden des Grundzustands ist ein NP-schweres Optimierungsproblem.

Voraussetzungen

Bevor du mit diesem Tutorial beginnst, stelle sicher, dass Folgendes installiert ist:

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

• Qiskit-Addon Optimization Mapper (pip install qiskit_addon_opt_mapper)

• NumPy (pip install numpy)

Du benötigst außerdem eine Berechtigung für den Zugriff auf die ParityQC Twine Optimizer-Funktion. Um Zugriff zu beantragen, fülle dieses Formular aus.

Setup

(Dieser Code setzt voraus, dass du dein Konto bereits in deiner lokalen Umgebung gespeichert hast.)

Importiere zunächst alle erforderlichen Pakete für dieses 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

Lade den Parity Twine Optimizer aus dem Qiskit Functions-Katalog:

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

Schritt 1: Das Problem als Zielfunktion definieren

Anstatt das SK-Problem aus einer Bibliothek zu beziehen, wie wir es beim Market-Split-Problem tun, formulierst du es direkt.

Die Funktion generate_sk_problem formuliert das SK-Problem direkt im erforderlichen Dictionary-Format. Die einzige erforderliche Eingabe ist n, die Anzahl der Spins im Modell.

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

Verwende die Funktion _generate_couplings_sk_model, um die zufälligen Kopplungsterme im SK-Modell zu berechnen. Für mehr Kontrolle über die Kopplungen kannst du optionale Argumente verwenden, die im Docstring der Funktion erklärt werden.

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

Schritt 2: Das Problem mit dem Parity Twine Optimizer lösen

Mit den obigen Funktionen kannst du das SK-Problem einrichten und mit dem Twine Optimizer und einem ausgewählten IBM Quantum®-Backend eine Lösung finden.

Um die Funktion auszuführen, wähle ein geeignetes Backend aus, zum Beispiel ibm_phoenix.

Optional kannst du über options die Übermittlung genauer steuern:

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

wobei shots eine Ganzzahl ist, die die Anzahl der Schaltkreisausführungen angibt, postprocessing_level bestimmt, ob eine Nachbearbeitung auf das Ergebnis angewendet wird, transpile_only legt fest, ob das Problem nur in einen Schaltkreis transpiliert (und nicht gelöst) wird, und job_tags ist die Bezeichnung, mit der der Job auf der IBM Quantum Platform identifiziert wird.

Die Größe des SK-Modells wird durch $N$ definiert, die Anzahl der wechselwirkenden Spins. Sobald du $N$ gewählt hast, generiert der obige Code das Problem für n_spins.

Führe den Optimizer aus:

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}")

Überprüfe den Job-Status:

# Monitor the job status
function_job.status()

Rufe die Ergebnisse ab:

result = function_job.result()

result

Das Ergebnis hat die Form:

{
    '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},
        },
    }
}

wobei das solution-Dictionary den im Problem definierten Qubits entspricht und deren optimierte Spin-Werte für den Hamiltonoperator des SK-Modells angibt. Diese spezifische Abfolge von Spins in der optimalen Lösung stellt die Konfiguration dar, die die Gesamtenergie des Systems basierend auf den gegebenen zufälligen Wechselwirkungsstärken minimiert. Im SK-Modell kann dies als der Zustand niedrigster Energie eines ungeordneten magnetischen Systems betrachtet werden.

metadata liefert Informationen zur Transpilation (Anzahl/Tiefe der Zwei-Qubit-Gates, verwendete Gates, aktive Qubits) und verschiedene Laufzeiten.

Nächste Schritte

Empfehlungen

• Fordere Zugriff auf die Funktion an, indem du dieses Formular ausfüllst.
• Besuche die API-Referenz für diese Qiskit-Funktion.
• Lies den Leitfaden.
• Probiere das Tutorial aus, um den Parity Twine Optimizer auf das Market-Split-Problem anzuwenden.
• Sieh dir den ArXiv-Preprint Connectivity-aware Synthesis of Quantum Algorithms, Drier et al. (2025) an.