Optimizare binară de ordin superior cu Optimization Solver de la Q-CTRL

Notă

Qiskit Functions este o funcționalitate experimentală disponibilă doar pentru utilizatorii planurilor IBM Quantum® Premium, Flex și On-Prem (prin API-ul IBM Quantum Platform). Aceasta se află în stadiu de previzualizare și poate fi modificată.

Estimare de utilizare: 24 de minute pe un procesor Heron r2. (NOTĂ: Aceasta este doar o estimare. Timpul tău de execuție poate varia.)

Fundal

Acest tutorial demonstrează cum să rezolvi o problemă de optimizare binară de ordin superior (HOBO) folosind Optimization Solver, o funcție Qiskit de la Q-CTRL Fire Opal. Exemplul prezentat în acest tutorial este o problemă de optimizare concepută pentru a găsi energia stării de bază a unui model Ising cu 156 qubiți și legături aleatoare, care conține termeni cubici. Optimization Solver poate fi utilizat pentru probleme de optimizare generale ce pot fi definite ca o funcție obiectiv.

Optimization Solver automatizează complet pașii de implementare conștienți de hardware pentru rezolvarea problemelor de optimizare pe hardware cuantic și, prin utilizarea Performance Management pentru execuția cuantică, obține soluții precise la scară utilă. Pentru o prezentare detaliată a fluxului de lucru complet al Optimization Solver și a rezultatelor de benchmarking, consultă manuscrisul publicat.

Acest tutorial parcurge următorii pași:

1. Definirea problemei ca funcție obiectiv
2. Rularea algoritmului hibrid folosind Fire Opal Optimization Solver
3. Evaluarea rezultatelor

Cerințe

Înainte de a începe acest tutorial, asigură-te că ai instalat următoarele:

• Qiskit Functions (pip install qiskit-ibm-catalog)
• SymPy (pip install sympy)

Va trebui, de asemenea, să obții acces la funcția Optimization Solver. Completează formularul pentru a solicita accesul.

Configurare

Mai întâi, importă pachetele și instrumentele necesare.

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

# Qiskit Functions Catalog
from qiskit_ibm_catalog import QiskitFunctionsCatalog

# SymPy tools for constructing objective function
from sympy import Poly
from sympy import symbols, srepr

# Tools for plotting and evaluating results
import numpy as np
import matplotlib.pyplot as plt
from sympy import lambdify

Definește-ți credențialele IBM Quantum Platform, care vor fi utilizate pe parcursul tutorialului pentru autentificarea la Qiskit Runtime și Qiskit Functions.

# Credentials
token = "<YOUR-API_KEY>"  # Use the 44-characters API_KEY you have created and saved from the IBM Quantum Platform Home dashboard
instance = "<YOUR_CRN>"

Pasul 1: Definirea problemei ca funcție obiectiv

Optimization Solver acceptă ca intrare fie o funcție obiectiv, fie un graf. În acest tutorial, problema de minimizare a sticlei de spin Ising este definită ca o funcție obiectiv și a fost adaptată pentru topologia heavy-hex a dispozitivelor IBM®.

Deoarece această funcție obiectiv conține termeni cubici, pătratici și liniari, se încadrează în clasa HOBO de probleme, cunoscute ca fiind considerabil mai dificile de rezolvat decât problemele convenționale de optimizare binară pătratică nerestricționată (QUBO).

Pentru o discuție detaliată despre construcția definiției problemei și rezultatele anterioare obținute cu Optimization Solver, consultă acest manuscris tehnic. Problema a fost inițial definită și evaluată în cadrul unui articol publicat de Los Alamos National Laboratory și a fost adaptată pentru a utiliza lățimea completă a procesoarelor IBM Quantum Heron cu 156 qubiți.

qubit_count = 156

# Create symbolic variables to represent qubits
x = symbols([f"x[{i}]" for i in range(qubit_count)])

# # Define a polynomial representing a spin glass model
spin_glass_poly = Poly(
    -4 * x[0] * x[1]
    - 8 * x[1] * x[2] * x[3]
    + 8 * x[1] * x[2]
    + 4 * x[1] * x[3]
    - 4 * x[2]
    + 8 * x[3] * x[4] * x[5]
    - 4 * x[3] * x[5]
    - 8 * x[3] * x[16] * x[23]
    + 4 * x[3] * x[23]
    - 2 * x[3]
    - 4 * x[4]
    - 8 * x[5] * x[6] * x[7]
    + 8 * x[5] * x[6]
    + 4 * x[5] * x[7]
    - 2 * x[5]
    + 8 * x[6] * x[7]
    - 4 * x[6]
    - 8 * x[7] * x[8] * x[9]
    + 4 * x[7] * x[9]
    - 8 * x[7] * x[17] * x[27]
    + 4 * x[7] * x[27]
    - 6 * x[7]
    + 8 * x[8] * x[9]
    + 8 * x[9] * x[10] * x[11]
    - 4 * x[9] * x[11]
    - 2 * x[9]
    - 8 * x[10] * x[11]
    + 4 * x[10]
    - 8 * x[11] * x[12] * x[13]
    + 4 * x[11] * x[13]
    - 8 * x[11] * x[18] * x[31]
    + 8 * x[11] * x[18]
    + 4 * x[11] * x[31]
    - 2 * x[11]
    + 8 * x[12] * x[13]
    + 8 * x[13] * x[14] * x[15]
    - 4 * x[13] * x[15]
    - 2 * x[13]
    - 8 * x[14] * x[15]
    + 4 * x[14]
    - 8 * x[15] * x[19] * x[35]
    + 8 * x[15] * x[19]
    + 4 * x[15] * x[35]
    - 2 * x[15]
    + 8 * x[16] * x[23]
    + 8 * x[17] * x[27]
    - 4 * x[17]
    + 8 * x[18] * x[31]
    - 8 * x[18]
    + 8 * x[19] * x[35]
    - 8 * x[19]
    + 4 * x[20] * x[21]
    - 4 * x[20]
    - 8 * x[21] * x[22] * x[23]
    + 8 * x[21] * x[22]
    + 4 * x[21] * x[23]
    - 8 * x[21] * x[36] * x[41]
    + 4 * x[21] * x[41]
    - 4 * x[21]
    + 8 * x[22] * x[23]
    - 8 * x[22]
    + 8 * x[23] * x[24] * x[25]
    - 4 * x[23] * x[25]
    - 10 * x[23]
    - 8 * x[24] * x[25]
    + 8 * x[25] * x[26] * x[27]
    - 8 * x[25] * x[26]
    - 4 * x[25] * x[27]
    + 8 * x[25] * x[37] * x[45]
    - 8 * x[25] * x[37]
    - 4 * x[25] * x[45]
    + 14 * x[25]
    - 8 * x[26] * x[27]
    + 4 * x[26]
    + 8 * x[27] * x[28] * x[29]
    - 4 * x[27] * x[29]
    - 2 * x[27]
    - 8 * x[28] * x[29]
    - 8 * x[29] * x[30] * x[31]
    + 4 * x[29] * x[31]
    + 8 * x[29] * x[38] * x[49]
    - 8 * x[29] * x[38]
    - 4 * x[29] * x[49]
    + 6 * x[29]
    + 8 * x[30] * x[31]
    - 4 * x[30]
    - 8 * x[31] * x[32] * x[33]
    + 4 * x[31] * x[33]
    - 6 * x[31]
    + 8 * x[33] * x[34] * x[35]
    - 4 * x[33] * x[35]
    - 8 * x[33] * x[39] * x[53]
    + 8 * x[33] * x[39]
    + 4 * x[33] * x[53]
    - 6 * x[33]
    - 8 * x[34] * x[35]
    + 2 * x[35]
    + 8 * x[36] * x[41]
    - 8 * x[37] * x[45]
    + 4 * x[37]
    - 8 * x[38] * x[49]
    + 4 * x[38]
    + 4 * x[40] * x[41]
    - 8 * x[41] * x[42] * x[43]
    + 4 * x[41] * x[43]
    - 8 * x[41]
    + 8 * x[42] * x[43]
    - 4 * x[42]
    - 8 * x[43] * x[44] * x[45]
    + 8 * x[43] * x[44]
    + 4 * x[43] * x[45]
    - 8 * x[43] * x[56] * x[63]
    + 4 * x[43] * x[63]
    - 6 * x[43]
    - 4 * x[44]
    - 8 * x[45] * x[46] * x[47]
    + 4 * x[45] * x[47]
    + 2 * x[45]
    + 4 * x[46]
    - 8 * x[47] * x[48] * x[49]
    + 8 * x[47] * x[48]
    + 4 * x[47] * x[49]
    - 8 * x[47] * x[57] * x[67]
    + 4 * x[47] * x[67]
    - 2 * x[47]
    - 4 * x[48]
    - 8 * x[49] * x[50] * x[51]
    + 8 * x[49] * x[50]
    + 4 * x[49] * x[51]
    - 2 * x[49]
    + 8 * x[50] * x[51]
    - 8 * x[50]
    - 8 * x[51] * x[52] * x[53]
    + 8 * x[51] * x[52]
    + 4 * x[51] * x[53]
    - 8 * x[51] * x[58] * x[71]
    + 4 * x[51] * x[71]
    - 6 * x[51]
    + 8 * x[52] * x[53]
    - 8 * x[52]
    + 8 * x[53] * x[54] * x[55]
    - 8 * x[53] * x[54]
    - 4 * x[53] * x[55]
    - 2 * x[53]
    + 4 * x[54]
    - 8 * x[55] * x[59] * x[75]
    + 4 * x[55] * x[75]
    - 2 * x[55]
    + 8 * x[56] * x[63]
    + 8 * x[57] * x[67]
    - 4 * x[57]
    + 8 * x[58] * x[71]
    + 8 * x[59] * x[75]
    - 4 * x[59]
    + 4 * x[60] * x[61]
    + 8 * x[61] * x[62] * x[63]
    - 4 * x[61] * x[63]
    + 8 * x[61] * x[76] * x[81]
    - 8 * x[61] * x[76]
    - 4 * x[61] * x[81]
    - 8 * x[63] * x[64] * x[65]
    + 8 * x[63] * x[64]
    + 4 * x[63] * x[65]
    - 6 * x[63]
    + 8 * x[65] * x[66] * x[67]
    - 8 * x[65] * x[66]
    - 4 * x[65] * x[67]
    - 8 * x[65] * x[77] * x[85]
    + 4 * x[65] * x[85]
    + 2 * x[65]
    + 4 * x[66]
    - 8 * x[67] * x[68] * x[69]
    + 8 * x[67] * x[68]
    + 4 * x[67] * x[69]
    - 10 * x[67]
    + 8 * x[68] * x[69]
    - 4 * x[68]
    + 8 * x[69] * x[70] * x[71]
    - 4 * x[69] * x[71]
    - 8 * x[69] * x[78] * x[89]
    + 4 * x[69] * x[89]
    - 6 * x[69]
    + 8 * x[71] * x[72] * x[73]
    - 8 * x[71] * x[72]
    - 4 * x[71] * x[73]
    + 2 * x[71]
    - 8 * x[72] * x[73]
    + 8 * x[72]
    - 8 * x[73] * x[74] * x[75]
    + 8 * x[73] * x[74]
    + 4 * x[73] * x[75]
    - 8 * x[73] * x[79] * x[93]
    + 8 * x[73] * x[79]
    + 4 * x[73] * x[93]
    - 6 * x[73]
    + 8 * x[74] * x[75]
    - 4 * x[74]
    - 10 * x[75]
    + 4 * x[76]
    + 8 * x[78] * x[89]
    - 4 * x[78]
    - 4 * x[79]
    - 4 * x[80] * x[81]
    + 4 * x[80]
    - 8 * x[81] * x[82] * x[83]
    + 8 * x[81] * x[82]
    + 4 * x[81] * x[83]
    + 8 * x[82] * x[83]
    - 8 * x[82]
    - 8 * x[83] * x[84] * x[85]
    + 4 * x[83] * x[85]
    - 8 * x[83] * x[96] * x[103]
    + 4 * x[83] * x[103]
    - 2 * x[83]
    - 8 * x[85] * x[86] * x[87]
    + 8 * x[85] * x[86]
    + 4 * x[85] * x[87]
    - 6 * x[85]
    + 8 * x[86] * x[87]
    - 4 * x[86]
    - 8 * x[87] * x[88] * x[89]
    + 4 * x[87] * x[89]
    + 8 * x[87] * x[97] * x[107]
    - 8 * x[87] * x[97]
    - 4 * x[87] * x[107]
    + 2 * x[87]
    + 4 * x[88]
    - 8 * x[89] * x[90] * x[91]
    + 8 * x[89] * x[90]
    + 4 * x[89] * x[91]
    - 10 * x[89]
    + 8 * x[90] * x[91]
    - 8 * x[90]
    - 8 * x[91] * x[92] * x[93]
    + 4 * x[91] * x[93]
    - 8 * x[91] * x[98] * x[111]
    + 8 * x[91] * x[98]
    + 4 * x[91] * x[111]
    - 10 * x[91]
    + 8 * x[92] * x[93]
    - 4 * x[92]
    - 8 * x[93] * x[94] * x[95]
    + 4 * x[93] * x[95]
    - 6 * x[93]
    + 8 * x[95] * x[99] * x[115]
    - 8 * x[95] * x[99]
    - 4 * x[95] * x[115]
    + 2 * x[95]
    + 4 * x[96]
    - 8 * x[97] * x[107]
    + 4 * x[97]
    - 4 * x[98]
    - 8 * x[99] * x[115]
    + 4 * x[99]
    - 4 * x[100] * x[101]
    + 8 * x[101] * x[102] * x[103]
    - 8 * x[101] * x[102]
    - 4 * x[101] * x[103]
    - 8 * x[101] * x[116] * x[121]
    + 8 * x[101] * x[116]
    + 4 * x[101] * x[121]
    + 4 * x[101]
    - 8 * x[103] * x[104] * x[105]
    + 4 * x[103] * x[105]
    + 2 * x[103]
    + 8 * x[105] * x[106] * x[107]
    - 4 * x[105] * x[107]
    - 8 * x[105] * x[117] * x[125]
    + 4 * x[105] * x[125]
    + 2 * x[105]
    - 8 * x[106] * x[107]
    + 4 * x[106]
    + 8 * x[107] * x[108] * x[109]
    - 4 * x[107] * x[109]
    + 6 * x[107]
    - 4 * x[108]
    + 8 * x[109] * x[110] * x[111]
    - 4 * x[109] * x[111]
    - 8 * x[109] * x[118] * x[129]
    + 4 * x[109] * x[129]
    + 2 * x[109]
    - 8 * x[110] * x[111]
    + 4 * x[110]
    - 8 * x[111] * x[112] * x[113]
    + 8 * x[111] * x[112]
    + 4 * x[111] * x[113]
    + 2 * x[111]
    + 8 * x[112] * x[113]
    - 8 * x[112]
    - 8 * x[113] * x[114] * x[115]
    + 4 * x[113] * x[115]
    - 8 * x[113] * x[119] * x[133]
    + 4 * x[113] * x[133]
    - 2 * x[113]
    + 6 * x[115]
    - 4 * x[116]
    + 4 * x[118]
    + 4 * x[119]
    + 4 * x[120] * x[121]
    - 8 * x[121] * x[122] * x[123]
    + 4 * x[121] * x[123]
    - 4 * x[121]
    + 4 * x[122]
    - 8 * x[123] * x[124] * x[125]
    + 4 * x[123] * x[125]
    - 8 * x[123] * x[136] * x[143]
    + 4 * x[123] * x[143]
    - 2 * x[123]
    + 8 * x[124] * x[125]
    - 4 * x[124]
    + 8 * x[125] * x[126] * x[127]
    - 8 * x[125] * x[126]
    - 4 * x[125] * x[127]
    + 2 * x[125]
    - 8 * x[127] * x[128] * x[129]
    + 8 * x[127] * x[128]
    + 4 * x[127] * x[129]
    + 8 * x[127] * x[137] * x[147]
    - 8 * x[127] * x[137]
    - 4 * x[127] * x[147]
    - 2 * x[127]
    + 8 * x[129] * x[130] * x[131]
    - 4 * x[129] * x[131]
    + 2 * x[129]
    - 4 * x[130]
    - 8 * x[131] * x[132] * x[133]
    + 4 * x[131] * x[133]
    - 8 * x[131] * x[138] * x[151]
    + 4 * x[131] * x[151]
    - 2 * x[131]
    + 8 * x[133] * x[134] * x[135]
    - 4 * x[133] * x[135]
    + 2 * x[133]
    - 8 * x[134] * x[135]
    + 4 * x[134]
    - 8 * x[135] * x[139] * x[155]
    + 8 * x[135] * x[139]
    + 4 * x[135] * x[155]
    + 2 * x[135]
    + 8 * x[136] * x[143]
    - 4 * x[136]
    + 4 * x[138]
    + 8 * x[139] * x[155]
    - 4 * x[139]
    - 4 * x[140] * x[141]
    - 8 * x[141] * x[142] * x[143]
    + 8 * x[141] * x[142]
    + 4 * x[141] * x[143]
    + 8 * x[142] * x[143]
    - 8 * x[142]
    - 8 * x[143] * x[144] * x[145]
    + 8 * x[143] * x[144]
    + 4 * x[143] * x[145]
    - 14 * x[143]
    + 8 * x[144] * x[145]
    - 8 * x[144]
    - 8 * x[145] * x[146] * x[147]
    + 8 * x[145] * x[146]
    + 4 * x[145] * x[147]
    - 6 * x[145]
    + 8 * x[146] * x[147]
    - 4 * x[146]
    - 8 * x[147] * x[148] * x[149]
    + 8 * x[147] * x[148]
    + 4 * x[147] * x[149]
    - 6 * x[147]
    - 4 * x[148]
    - 8 * x[149] * x[150] * x[151]
    + 8 * x[149] * x[150]
    + 4 * x[149] * x[151]
    - 6 * x[149]
    + 8 * x[151] * x[152] * x[153]
    - 4 * x[151] * x[153]
    + 2 * x[151]
    + 8 * x[153] * x[154] * x[155]
    - 8 * x[153] * x[154]
    - 4 * x[153] * x[155]
    + 2 * x[153]
    - 8 * x[154] * x[155]
    + 4 * x[154]
    - 2 * x[155]
    + 46,
    x,
    domain="ZZ",
)

Pasul 2: Rulează algoritmul hibrid folosind Fire Opal Optimization Solver

Acum folosește Qiskit Function Optimization Solver pentru a rula algoritmul. În fundal, Optimization Solver se ocupă de maparea problemei la un algoritm cuantic hibrid, rularea circuitelor cuantice cu supresie a erorilor și efectuarea optimizării clasice.

# Authenticate to the Qiskit Functions Catalog
catalog = QiskitFunctionsCatalog(
    token=token,
    instance=instance,
)

# Load the function
solver = catalog.load("q-ctrl/optimization_solver")

Verifică că device-ul ales are cel puțin 156 de qubiți.

# Specify the target backend name
backend_name = "<CHOOSE_A_BACKEND>"

Solver-ul acceptă o reprezentare sub formă de șir de caractere a funcției obiectiv.

# Convert the objective function to string format
spin_glass_poly_as_str = srepr(spin_glass_poly)

# Run the problem
spin_glass_job = solver.run(
    problem=spin_glass_poly_as_str,
    run_options={"backend_name": backend_name},
)

Poți folosi familiarele API-uri Qiskit Serverless pentru a verifica statusul sarcinii Qiskit Function:

# Get job status
spin_glass_job.status()

Solver-ul returnează un dicționar cu soluția și metadatele asociate, cum ar fi șirul de biți al soluției, numărul de iterații și maparea variabilelor la șirul de biți. Pentru o definiție completă a intrărilor și ieșirilor Solver-ului, consultă documentația.

# Poll for results
result = spin_glass_job.result()

# Get the final bitstring distribution and set the number of shots
distribution = result["final_bitstring_distribution"]

Pasul 3: Evaluează rezultatele

# Get the solution ground state energy
print(f"Minimum ground state energy: {result["solution_bitstring_cost"]}")

Minimum ground state energy: -242.0

Solver-ul a găsit soluția corectă, care a fost validată folosind software de optimizare clasic. Complexitatea acestei probleme la scară de utilitate necesită un software de optimizare avansat pentru a fi rezolvată clasic, cum ar fi IBM ILOG CPLEX Optimization Studio (CPLEX) sau Gurobi Optimization. Ca analiză vizuală a calității rezultatelor, poți reprezenta grafic rezultatele calculând valorile costurilor din șirurile de biți și probabilitățile acestora. Pentru comparație, trasează rezultatele alături de o distribuție de șiruri de biți eșantionate aleatoriu, ceea ce este echivalent cu o soluție clasică de tip „forță brută". Dacă algoritmul găsește în mod constant costuri mai mici, aceasta sugerează că algoritmul cuantic rezolvă eficient problema de optimizare.

def plot_cost_histogram(
    costs, probabilities, distribution, qubit_count, bitstring_cost
):
    """Plots a histogram comparing the cost distributions of Q-CTRL Solver and random sampling."""

    # Set figure DPI for higher resolution and font size for labels
    plt.rcParams["figure.dpi"] = 300
    plt.rcParams.update({"font.size": 6})  # Set default font size to 6

    # Define labels and colors for the plot
    labels = ["Q-CTRL Solver", "Random Sampling"]
    colors = ["#680CE9", "#E04542"]

    # Calculate total shots (total number of bitstrings in the distribution)
    shots = sum(distribution.values())

    # Generate random bitstrings for comparison (random sampling)
    rng = np.random.default_rng(seed=0)
    random_array = rng.integers(
        0, 2, size=(shots, qubit_count)
    )  # Generate random bitstrings (0 or 1 for each qubit)
    random_bitstrings = ["".join(row.astype(str)) for row in random_array]

    # Compute the cost for each random bitstring
    random_costs = [bitstring_cost(k) for k in random_bitstrings]

    # Set uniform probabilities for the random sampling
    random_probabilities = (
        np.ones(shape=(shots,)) / shots
    )  # Equal probability for each random bitstring

    # Find the minimum and maximum costs for binning the histogram
    min_cost = np.min(costs)
    max_cost = np.max(random_costs)

    # Create a histogram plot with a smaller figure size (4x2 inches)
    fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(4, 2))

    # Plot histograms for the Q-CTRL solver and random sampling costs
    _, _, _ = ax.hist(
        [costs, random_costs],  # Data for the two histograms
        np.arange(min_cost, max_cost, 2),  # Bins for the histogram
        weights=[
            probabilities,
            random_probabilities,
        ],  # Probabilities for each data set
        label=labels,  # Labels for the legend
        color=colors,  # Colors for each histogram
        histtype="stepfilled",  # Filled step histogram
        align="mid",  # Align bars to the bin center
        alpha=0.8,  # Transparency
    )

    # Set the x and y labels for the plot
    ax.set_xlabel("Cost")
    ax.set_ylabel("Probability")

    # Add the legend to the plot
    ax.legend()

    # Show the plot
    plt.show()

# Convert spin_glass_poly into a NumPy-compatible function
poly_as_numpy_function = lambdify(x, spin_glass_poly.as_expr(), "numpy")

# Function to compute the cost of a given bitstring using spin_glass_poly
def bitstring_cost(bitstring: str) -> float:
    # Convert bitstring to a reversed list of integers (0s and 1s)
    return float(
        poly_as_numpy_function(*[int(b) for b in str(bitstring[::-1])])
    )

# Calculate the cost of each bitstring in the distribution
costs = [bitstring_cost(k) for k, _ in distribution.items()]

# Extract probabilities from the bitstring distribution
probabilities = np.array([v for _, v in distribution.items()])
probabilities = probabilities / sum(
    probabilities
)  # Normalize to get probabilities

plot_cost_histogram(
    costs, probabilities, distribution, qubit_count, bitstring_cost
)

Deoarece scopul acestui algoritm de optimizare este de a găsi starea fundamentală minimă a modelului Ising, valorile mai mici indică soluții mai bune. Prin urmare, este vizibil că soluțiile generate de Fire Opal Optimization Solver depășesc cu mult selecția aleatorie.