Setup¶
Google Colab¶
Click the badge above to open this notebook in Colab. The notebook installs or activates dependencies in the setup cells below.
Local installation¶
Run the following from the repository root before opening this notebook locally:
python -m pip install dimod dwave-neal matplotlib networkx numpy pandas pyomo scipyGLPK is used for the integer-programming comparison cells.
GLPK (required for ILP sections):
macOS:
brew install glpkUbuntu/Debian:
sudo apt-get install glpk-utilsWindows: download from http://
winglpk .sourceforge .net/
Learning objectives¶
By the end of this notebook you will be able to:
Describe QUBO and Ising model forms and the role of linear, quadratic, and offset terms.
Convert constrained binary optimization problems into QUBO models using penalty terms.
Build Binary Quadratic Models with
dimodand solve them with simulated annealing.Validate sampled solutions and connect QUBO penalties to graph-coloring feasibility.
Prerequisites¶
Mathematical background: Binary variables, matrix notation, quadratic objectives, and basic graph terminology.
Prior notebooks: Notebook 1 (MathProg) or equivalent experience with binary integer programming.
Accounts required: None; all examples use local Python packages.
Python version: Python 3.9+ with the D-Wave Ocean packages used by this repository.
Quadratic Unconstrained Binary Optimization¶
This notebook explains the basics of QUBO modeling. We use D-Wave’s dimod package to build QUBOs and neal to solve them with simulated annealing. We also use SymPy for the symbolic computation behind the Gröbner-basis examples and NetworkX to represent graph problems.
QUBO problem statement¶
We define a QUBO as the following optimization problem:
where we optimize over binary variables , on a constrained graph defined by an adjacency matrix . We also include an arbitrary offset .
QUBO example¶
Suppose we want to solve the following problem via QUBO $$ \min_{\mathbf{x}} 2𝑥_0+4𝑥_1+4𝑥_2+4𝑥_3+4𝑥_4+4𝑥_5+5𝑥_6+4𝑥_7+5𝑥_8+6𝑥_9+5𝑥_{10} \ s.t. \begin{bmatrix} 1 & 0 & 0 & 1 & 1 & 1 & 0 & 1 & 1 & 1 & 1\ 0 & 1 & 0 & 1 & 0 & 1 & 1 & 0 & 1 & 1 & 1\ 0 & 0 & 1 & 0 & 1 & 0 & 1 & 1 & 1 & 1 & 1 \end{bmatrix}\mathbf{x}=
\mathbf{x} \in {0,1 }^{11} $$
Notebook Cell
# If using this on Google Colab, we need to install the packages
try:
import google.colab
IN_COLAB = True
except:
IN_COLAB = False
# Install the packages used by this notebook
if IN_COLAB:
!pip install -q pyomo dimod dwave-neal scipy pandas networkx matplotlib# Import the Pyomo library, which can be installed via pip, conda, or from GitHub: https://github.com/Pyomo/pyomo
import pyomo.environ as pyo
# Import the Dwave packages dimod and neal
import dimod
import neal
# Import Matplotlib to generate plots
import matplotlib.pyplot as plt
# Import numpy and scipy for certain numerical calculations below
import numpy as np
from collections import Counter
import pandas as pd
import networkx as nxFirst we would write this problem as an unconstrained one by penalizing the linear constraints as quadratics in the objective. Let’s first define the problem parameters
A = np.array([[1, 0, 0, 1, 1, 1, 0, 1, 1, 1, 1],
[0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1],
[0, 0, 1, 0, 1, 0, 1, 1, 1, 1, 1]])
b = np.array([1, 1, 1])
c = np.array([2, 4, 4, 4, 4, 4, 5, 4, 5,6, 5])
In order to define the matrix, we first write the problem
as follows:
Exploiting the fact that for , we can make the linear terms appear in the diagonal of the matrix.
Penalty parameter rationale¶
The penalty term must be large enough that any infeasible assignment is worse than the objective improvement it might gain by violating the constraint. For a binary constraint A x = b, the residual A x - b is integer-valued; the smallest nonzero violation has squared penalty at least 1. A conservative sufficient bound is rho = sum(abs(c)) + epsilon, because changing binary variables can improve the linear objective by at most sum(abs(c_i)).
Worked example: if sum(abs(c_i)) = 6 and rho = 5.9, an infeasible assignment with violation 1 can gain 6 objective units while paying only 5.9 penalty units, so it can look 0.1 units better than a feasible assignment. Choosing rho = 6 + epsilon closes that gap for this bounded binary model.
This is a sufficient bound for this model family, not a universal rule. Constraints with non-binary variables, non-integer residuals, or a larger objective range need a problem-specific penalty analysis. See Glover, Kochenberger, and Du (2019), “A Tutorial on Formulating and Using QUBO Models.”
epsilon = 1
rho = np.sum(np.abs(c)) + epsilon
Q = rho*np.matmul(A.T,A)
Q += np.diag(c)
Q -= rho*2*np.diag(np.matmul(b.T,A))
Beta = rho*np.matmul(b.T,b)
print(Q)
print(Beta)[[ -46 0 0 48 48 48 0 48 48 48 48]
[ 0 -44 0 48 0 48 48 0 48 48 48]
[ 0 0 -44 0 48 0 48 48 48 48 48]
[ 48 48 0 -92 48 96 48 48 96 96 96]
[ 48 0 48 48 -92 48 48 96 96 96 96]
[ 48 48 0 96 48 -92 48 48 96 96 96]
[ 0 48 48 48 48 48 -91 48 96 96 96]
[ 48 0 48 48 96 48 48 -92 96 96 96]
[ 48 48 48 96 96 96 96 96 -139 144 144]
[ 48 48 48 96 96 96 96 96 144 -138 144]
[ 48 48 48 96 96 96 96 96 144 144 -139]]
144
We can visualize the graph that defines this instance using the Q matrix as the adjacency matrix of a graph.
G = nx.from_numpy_array(Q)
# A fixed layout seed makes the drawing repeatable.
nx.draw(G, pos=nx.spring_layout(G, seed=314159), with_labels=True)
Let’s define a QUBO model and then solve it using D-Wave’s code for complete enumeration and simulated annealing (eventually with quantum annealing too!).
model = dimod.BinaryQuadraticModel.from_qubo(Q, offset=Beta)
print(model)BinaryQuadraticModel({0: -46.0, 1: -44.0, 2: -44.0, 3: -92.0, 4: -92.0, 5: -92.0, 6: -91.0, 7: -92.0, 8: -139.0, 9: -138.0, 10: -139.0}, {(3, 0): 96.0, (3, 1): 96.0, (4, 0): 96.0, (4, 2): 96.0, (4, 3): 96.0, (5, 0): 96.0, (5, 1): 96.0, (5, 3): 192.0, (5, 4): 96.0, (6, 1): 96.0, (6, 2): 96.0, (6, 3): 96.0, (6, 4): 96.0, (6, 5): 96.0, (7, 0): 96.0, (7, 2): 96.0, (7, 3): 96.0, (7, 4): 192.0, (7, 5): 96.0, (7, 6): 96.0, (8, 0): 96.0, (8, 1): 96.0, (8, 2): 96.0, (8, 3): 192.0, (8, 4): 192.0, (8, 5): 192.0, (8, 6): 192.0, (8, 7): 192.0, (9, 0): 96.0, (9, 1): 96.0, (9, 2): 96.0, (9, 3): 192.0, (9, 4): 192.0, (9, 5): 192.0, (9, 6): 192.0, (9, 7): 192.0, (9, 8): 288.0, (10, 0): 96.0, (10, 1): 96.0, (10, 2): 96.0, (10, 3): 192.0, (10, 4): 192.0, (10, 5): 192.0, (10, 6): 192.0, (10, 7): 192.0, (10, 8): 288.0, (10, 9): 288.0}, 144.0, 'BINARY')
Since the problem is relatively small (11 variables, combinations), we can afford to enumerate all the solutions.
exactSampler = dimod.reference.samplers.ExactSolver()
exactSamples = exactSampler.sample(model)# Some useful functions to get plots
def sparse_tick_positions(count, max_ticks=10):
if count <= max_ticks:
return np.arange(count, dtype=int)
return np.unique(np.linspace(0, count - 1, num=max_ticks, dtype=int))
def binary_sample_label(sample, variables, max_length=18):
bitstring = ''.join(str(int(sample[variable])) for variable in variables)
if len(bitstring) <= max_length:
return bitstring
side_length = (max_length - 3) // 2
return f'{bitstring[:side_length]}...{bitstring[-side_length:]}'
def plot_enumerate(results, title=None, max_xticks=10):
records = list(results.data(['sample', 'energy'], sorted_by='energy'))
energies = [datum.energy for datum in records]
positions = np.arange(len(records))
tick_positions = sparse_tick_positions(len(records), max_xticks)
fig, ax = plt.subplots(figsize=(10, 4))
ax.bar(positions, energies)
if results.vartype == dimod.BINARY:
variables = list(results.variables)
sample_labels = [binary_sample_label(datum.sample, variables) for datum in records]
tick_labels = [sample_labels[index] for index in tick_positions]
ax.set_xlabel('solution rank (selected bitstrings)')
ax.set_xticks(tick_positions, tick_labels, rotation=45, ha='right')
else:
ax.set_xlabel('solution rank')
ax.set_xticks(tick_positions, tick_positions)
ax.set_ylabel('Energy')
ax.set_title(str(title))
fig.tight_layout()
plt.show()
print("minimum energy:", min(energies))
def plot_samples(results, title=None):
plt.figure()
energies = [datum.energy for datum in results.data(
['energy'], sorted_by='energy')]
if results.vartype == dimod.BINARY:
samples = [''.join(c for c in str(datum.sample.values()).strip(
', ') if c.isdigit()) for datum in results.data(['sample'], sorted_by=None)]
plt.xlabel('bitstring for solution')
else:
samples = np.arange(len(energies))
plt.xlabel('solution')
counts = Counter(samples)
total = len(samples)
for key in counts:
counts[key] /= total
df = pd.DataFrame.from_dict(counts, orient='index').sort_index()
df.plot(kind='bar', legend=None)
plt.xticks(rotation=80)
plt.ylabel('Probabilities')
plt.title(str(title))
plt.show()
print("minimum energy:", min(energies))
def plot_energies(results, title=None, max_xticks=10):
counts = Counter()
for datum in results.data(['energy', 'num_occurrences'], sorted_by='energy'):
counts[datum.energy] += datum.num_occurrences
total = sum(counts.values())
energy_values = sorted(counts)
probabilities = [counts[energy] / total for energy in energy_values]
positions = np.arange(len(energy_values))
tick_positions = sparse_tick_positions(len(energy_values), max_xticks)
tick_labels = [f'{energy_values[index]:g}' for index in tick_positions]
fig, ax = plt.subplots(figsize=(10, 4))
ax.bar(positions, probabilities)
ax.set_xticks(tick_positions, tick_labels, rotation=45, ha='right')
ax.set_xlabel('Energy')
ax.set_ylabel('Probabilities')
ax.set_title(str(title))
fig.tight_layout()
plt.show()
print("minimum energy:", min(energy_values))plot_enumerate(exactSamples, title='Enumerate all solutions')
plot_energies(exactSamples, title='Enumerate all solutions')
minimum energy: 5.0
minimum energy: 5.0
Let’s now solve this QUBO via traditional Integer Programming.
# We do not need to worry about the tranformation to QUBO since dimod takes care of it
Q, c = model.to_qubo()
# Define the model
model_pyo = pyo.ConcreteModel(name='QUBO example as an IP, 47-779/785 QuIPML')
I = range(len(model))
J_idx = range(len(model))
#Define the original variables
model_pyo.x = pyo.Var(I, domain=pyo.Binary)
# Define the edges variables
model_pyo.y = pyo.Var(I, J_idx, domain=pyo.Binary)
obj_expr = c
# add model constraints
model_pyo.c1 = pyo.ConstraintList()
model_pyo.c2 = pyo.ConstraintList()
model_pyo.c3 = pyo.ConstraintList()
for (i,j) in Q.keys():
if i != j:
model_pyo.c1.add(model_pyo.y[i,j] >= model_pyo.x[i] + model_pyo.x[j] - 1)
model_pyo.c2.add(model_pyo.y[i,j] <= model_pyo.x[i])
model_pyo.c3.add(model_pyo.y[i,j] <= model_pyo.x[j])
obj_expr += Q[i,j]*model_pyo.y[i,j]
else:
obj_expr += Q[i,j]*model_pyo.x[i]
# Define the objective function
model_pyo.objective = pyo.Objective(expr = obj_expr, sense=pyo.minimize)
# Print the model
model_pyo.display()Model 'QUBO example as an IP, 47-779/785 QuIPML'
Variables:
x : Size=11, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
0 : 0 : None : 1 : False : True : Binary
1 : 0 : None : 1 : False : True : Binary
2 : 0 : None : 1 : False : True : Binary
3 : 0 : None : 1 : False : True : Binary
4 : 0 : None : 1 : False : True : Binary
5 : 0 : None : 1 : False : True : Binary
6 : 0 : None : 1 : False : True : Binary
7 : 0 : None : 1 : False : True : Binary
8 : 0 : None : 1 : False : True : Binary
9 : 0 : None : 1 : False : True : Binary
10 : 0 : None : 1 : False : True : Binary
y : Size=121, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}*{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
(0, 0) : 0 : None : 1 : False : True : Binary
(0, 1) : 0 : None : 1 : False : True : Binary
(0, 2) : 0 : None : 1 : False : True : Binary
(0, 3) : 0 : None : 1 : False : True : Binary
(0, 4) : 0 : None : 1 : False : True : Binary
(0, 5) : 0 : None : 1 : False : True : Binary
(0, 6) : 0 : None : 1 : False : True : Binary
(0, 7) : 0 : None : 1 : False : True : Binary
(0, 8) : 0 : None : 1 : False : True : Binary
(0, 9) : 0 : None : 1 : False : True : Binary
(0, 10) : 0 : None : 1 : False : True : Binary
(1, 0) : 0 : None : 1 : False : True : Binary
(1, 1) : 0 : None : 1 : False : True : Binary
(1, 2) : 0 : None : 1 : False : True : Binary
(1, 3) : 0 : None : 1 : False : True : Binary
(1, 4) : 0 : None : 1 : False : True : Binary
(1, 5) : 0 : None : 1 : False : True : Binary
(1, 6) : 0 : None : 1 : False : True : Binary
(1, 7) : 0 : None : 1 : False : True : Binary
(1, 8) : 0 : None : 1 : False : True : Binary
(1, 9) : 0 : None : 1 : False : True : Binary
(1, 10) : 0 : None : 1 : False : True : Binary
(2, 0) : 0 : None : 1 : False : True : Binary
(2, 1) : 0 : None : 1 : False : True : Binary
(2, 2) : 0 : None : 1 : False : True : Binary
(2, 3) : 0 : None : 1 : False : True : Binary
(2, 4) : 0 : None : 1 : False : True : Binary
(2, 5) : 0 : None : 1 : False : True : Binary
(2, 6) : 0 : None : 1 : False : True : Binary
(2, 7) : 0 : None : 1 : False : True : Binary
(2, 8) : 0 : None : 1 : False : True : Binary
(2, 9) : 0 : None : 1 : False : True : Binary
(2, 10) : 0 : None : 1 : False : True : Binary
(3, 0) : 0 : None : 1 : False : True : Binary
(3, 1) : 0 : None : 1 : False : True : Binary
(3, 2) : 0 : None : 1 : False : True : Binary
(3, 3) : 0 : None : 1 : False : True : Binary
(3, 4) : 0 : None : 1 : False : True : Binary
(3, 5) : 0 : None : 1 : False : True : Binary
(3, 6) : 0 : None : 1 : False : True : Binary
(3, 7) : 0 : None : 1 : False : True : Binary
(3, 8) : 0 : None : 1 : False : True : Binary
(3, 9) : 0 : None : 1 : False : True : Binary
(3, 10) : 0 : None : 1 : False : True : Binary
(4, 0) : 0 : None : 1 : False : True : Binary
(4, 1) : 0 : None : 1 : False : True : Binary
(4, 2) : 0 : None : 1 : False : True : Binary
(4, 3) : 0 : None : 1 : False : True : Binary
(4, 4) : 0 : None : 1 : False : True : Binary
(4, 5) : 0 : None : 1 : False : True : Binary
(4, 6) : 0 : None : 1 : False : True : Binary
(4, 7) : 0 : None : 1 : False : True : Binary
(4, 8) : 0 : None : 1 : False : True : Binary
(4, 9) : 0 : None : 1 : False : True : Binary
(4, 10) : 0 : None : 1 : False : True : Binary
(5, 0) : 0 : None : 1 : False : True : Binary
(5, 1) : 0 : None : 1 : False : True : Binary
(5, 2) : 0 : None : 1 : False : True : Binary
(5, 3) : 0 : None : 1 : False : True : Binary
(5, 4) : 0 : None : 1 : False : True : Binary
(5, 5) : 0 : None : 1 : False : True : Binary
(5, 6) : 0 : None : 1 : False : True : Binary
(5, 7) : 0 : None : 1 : False : True : Binary
(5, 8) : 0 : None : 1 : False : True : Binary
(5, 9) : 0 : None : 1 : False : True : Binary
(5, 10) : 0 : None : 1 : False : True : Binary
(6, 0) : 0 : None : 1 : False : True : Binary
(6, 1) : 0 : None : 1 : False : True : Binary
(6, 2) : 0 : None : 1 : False : True : Binary
(6, 3) : 0 : None : 1 : False : True : Binary
(6, 4) : 0 : None : 1 : False : True : Binary
(6, 5) : 0 : None : 1 : False : True : Binary
(6, 6) : 0 : None : 1 : False : True : Binary
(6, 7) : 0 : None : 1 : False : True : Binary
(6, 8) : 0 : None : 1 : False : True : Binary
(6, 9) : 0 : None : 1 : False : True : Binary
(6, 10) : 0 : None : 1 : False : True : Binary
(7, 0) : 0 : None : 1 : False : True : Binary
(7, 1) : 0 : None : 1 : False : True : Binary
(7, 2) : 0 : None : 1 : False : True : Binary
(7, 3) : 0 : None : 1 : False : True : Binary
(7, 4) : 0 : None : 1 : False : True : Binary
(7, 5) : 0 : None : 1 : False : True : Binary
(7, 6) : 0 : None : 1 : False : True : Binary
(7, 7) : 0 : None : 1 : False : True : Binary
(7, 8) : 0 : None : 1 : False : True : Binary
(7, 9) : 0 : None : 1 : False : True : Binary
(7, 10) : 0 : None : 1 : False : True : Binary
(8, 0) : 0 : None : 1 : False : True : Binary
(8, 1) : 0 : None : 1 : False : True : Binary
(8, 2) : 0 : None : 1 : False : True : Binary
(8, 3) : 0 : None : 1 : False : True : Binary
(8, 4) : 0 : None : 1 : False : True : Binary
(8, 5) : 0 : None : 1 : False : True : Binary
(8, 6) : 0 : None : 1 : False : True : Binary
(8, 7) : 0 : None : 1 : False : True : Binary
(8, 8) : 0 : None : 1 : False : True : Binary
(8, 9) : 0 : None : 1 : False : True : Binary
(8, 10) : 0 : None : 1 : False : True : Binary
(9, 0) : 0 : None : 1 : False : True : Binary
(9, 1) : 0 : None : 1 : False : True : Binary
(9, 2) : 0 : None : 1 : False : True : Binary
(9, 3) : 0 : None : 1 : False : True : Binary
(9, 4) : 0 : None : 1 : False : True : Binary
(9, 5) : 0 : None : 1 : False : True : Binary
(9, 6) : 0 : None : 1 : False : True : Binary
(9, 7) : 0 : None : 1 : False : True : Binary
(9, 8) : 0 : None : 1 : False : True : Binary
(9, 9) : 0 : None : 1 : False : True : Binary
(9, 10) : 0 : None : 1 : False : True : Binary
(10, 0) : 0 : None : 1 : False : True : Binary
(10, 1) : 0 : None : 1 : False : True : Binary
(10, 2) : 0 : None : 1 : False : True : Binary
(10, 3) : 0 : None : 1 : False : True : Binary
(10, 4) : 0 : None : 1 : False : True : Binary
(10, 5) : 0 : None : 1 : False : True : Binary
(10, 6) : 0 : None : 1 : False : True : Binary
(10, 7) : 0 : None : 1 : False : True : Binary
(10, 8) : 0 : None : 1 : False : True : Binary
(10, 9) : 0 : None : 1 : False : True : Binary
(10, 10) : 0 : None : 1 : False : True : Binary
Objectives:
objective : Size=1, Index=None, Active=True
ERROR: evaluating object as numeric value: y[3,0]
(object: <class 'pyomo.core.base.var.VarData'>)
No value for uninitialized VarData object y[3,0]
ERROR: evaluating object as numeric value: objective
(object: <class 'pyomo.core.base.objective.ScalarObjective'>)
No value for uninitialized VarData object y[3,0]
Key : Active : Value
None : None : None
Constraints:
c1 : Size=47
Key : Lower : Body : Upper
1 : None : None : 0.0
2 : None : None : 0.0
3 : None : None : 0.0
4 : None : None : 0.0
5 : None : None : 0.0
6 : None : None : 0.0
7 : None : None : 0.0
8 : None : None : 0.0
9 : None : None : 0.0
10 : None : None : 0.0
11 : None : None : 0.0
12 : None : None : 0.0
13 : None : None : 0.0
14 : None : None : 0.0
15 : None : None : 0.0
16 : None : None : 0.0
17 : None : None : 0.0
18 : None : None : 0.0
19 : None : None : 0.0
20 : None : None : 0.0
21 : None : None : 0.0
22 : None : None : 0.0
23 : None : None : 0.0
24 : None : None : 0.0
25 : None : None : 0.0
26 : None : None : 0.0
27 : None : None : 0.0
28 : None : None : 0.0
29 : None : None : 0.0
30 : None : None : 0.0
31 : None : None : 0.0
32 : None : None : 0.0
33 : None : None : 0.0
34 : None : None : 0.0
35 : None : None : 0.0
36 : None : None : 0.0
37 : None : None : 0.0
38 : None : None : 0.0
39 : None : None : 0.0
40 : None : None : 0.0
41 : None : None : 0.0
42 : None : None : 0.0
43 : None : None : 0.0
44 : None : None : 0.0
45 : None : None : 0.0
46 : None : None : 0.0
47 : None : None : 0.0
c2 : Size=47
Key : Lower : Body : Upper
1 : None : None : 0.0
2 : None : None : 0.0
3 : None : None : 0.0
4 : None : None : 0.0
5 : None : None : 0.0
6 : None : None : 0.0
7 : None : None : 0.0
8 : None : None : 0.0
9 : None : None : 0.0
10 : None : None : 0.0
11 : None : None : 0.0
12 : None : None : 0.0
13 : None : None : 0.0
14 : None : None : 0.0
15 : None : None : 0.0
16 : None : None : 0.0
17 : None : None : 0.0
18 : None : None : 0.0
19 : None : None : 0.0
20 : None : None : 0.0
21 : None : None : 0.0
22 : None : None : 0.0
23 : None : None : 0.0
24 : None : None : 0.0
25 : None : None : 0.0
26 : None : None : 0.0
27 : None : None : 0.0
28 : None : None : 0.0
29 : None : None : 0.0
30 : None : None : 0.0
31 : None : None : 0.0
32 : None : None : 0.0
33 : None : None : 0.0
34 : None : None : 0.0
35 : None : None : 0.0
36 : None : None : 0.0
37 : None : None : 0.0
38 : None : None : 0.0
39 : None : None : 0.0
40 : None : None : 0.0
41 : None : None : 0.0
42 : None : None : 0.0
43 : None : None : 0.0
44 : None : None : 0.0
45 : None : None : 0.0
46 : None : None : 0.0
47 : None : None : 0.0
c3 : Size=47
Key : Lower : Body : Upper
1 : None : None : 0.0
2 : None : None : 0.0
3 : None : None : 0.0
4 : None : None : 0.0
5 : None : None : 0.0
6 : None : None : 0.0
7 : None : None : 0.0
8 : None : None : 0.0
9 : None : None : 0.0
10 : None : None : 0.0
11 : None : None : 0.0
12 : None : None : 0.0
13 : None : None : 0.0
14 : None : None : 0.0
15 : None : None : 0.0
16 : None : None : 0.0
17 : None : None : 0.0
18 : None : None : 0.0
19 : None : None : 0.0
20 : None : None : 0.0
21 : None : None : 0.0
22 : None : None : 0.0
23 : None : None : 0.0
24 : None : None : 0.0
25 : None : None : 0.0
26 : None : None : 0.0
27 : None : None : 0.0
28 : None : None : 0.0
29 : None : None : 0.0
30 : None : None : 0.0
31 : None : None : 0.0
32 : None : None : 0.0
33 : None : None : 0.0
34 : None : None : 0.0
35 : None : None : 0.0
36 : None : None : 0.0
37 : None : None : 0.0
38 : None : None : 0.0
39 : None : None : 0.0
40 : None : None : 0.0
41 : None : None : 0.0
42 : None : None : 0.0
43 : None : None : 0.0
44 : None : None : 0.0
45 : None : None : 0.0
46 : None : None : 0.0
47 : None : None : 0.0
Let’s install the MIP solver GLPK
Notebook Cell
# Let's install the LP/MIP solver GLPK
if IN_COLAB:
!apt-get install -y -qq glpk-utils# Define the solver GLPK
if IN_COLAB:
opt_glpk = pyo.SolverFactory('glpk', executable='/usr/bin/glpsol')
else:
opt_glpk = pyo.SolverFactory('glpk')
# Here we could use another solver, e.g. gurobi or cplex
# opt_gurobi = pyo.SolverFactory('gurobi')# We obtain the solution with GLPK
result_obj = opt_glpk.solve(model_pyo, tee=False)
model_pyo.display()Model 'QUBO example as an IP, 47-779/785 QuIPML'
Variables:
x : Size=11, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
0 : 0 : 0.0 : 1 : False : False : Binary
1 : 0 : 0.0 : 1 : False : False : Binary
2 : 0 : 0.0 : 1 : False : False : Binary
3 : 0 : 0.0 : 1 : False : False : Binary
4 : 0 : 0.0 : 1 : False : False : Binary
5 : 0 : 0.0 : 1 : False : False : Binary
6 : 0 : 0.0 : 1 : False : False : Binary
7 : 0 : 0.0 : 1 : False : False : Binary
8 : 0 : 0.0 : 1 : False : False : Binary
9 : 0 : 0.0 : 1 : False : False : Binary
10 : 0 : 1.0 : 1 : False : False : Binary
y : Size=121, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}*{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
(0, 0) : 0 : None : 1 : False : True : Binary
(0, 1) : 0 : None : 1 : False : True : Binary
(0, 2) : 0 : None : 1 : False : True : Binary
(0, 3) : 0 : None : 1 : False : True : Binary
(0, 4) : 0 : None : 1 : False : True : Binary
(0, 5) : 0 : None : 1 : False : True : Binary
(0, 6) : 0 : None : 1 : False : True : Binary
(0, 7) : 0 : None : 1 : False : True : Binary
(0, 8) : 0 : None : 1 : False : True : Binary
(0, 9) : 0 : None : 1 : False : True : Binary
(0, 10) : 0 : None : 1 : False : True : Binary
(1, 0) : 0 : None : 1 : False : True : Binary
(1, 1) : 0 : None : 1 : False : True : Binary
(1, 2) : 0 : None : 1 : False : True : Binary
(1, 3) : 0 : None : 1 : False : True : Binary
(1, 4) : 0 : None : 1 : False : True : Binary
(1, 5) : 0 : None : 1 : False : True : Binary
(1, 6) : 0 : None : 1 : False : True : Binary
(1, 7) : 0 : None : 1 : False : True : Binary
(1, 8) : 0 : None : 1 : False : True : Binary
(1, 9) : 0 : None : 1 : False : True : Binary
(1, 10) : 0 : None : 1 : False : True : Binary
(2, 0) : 0 : None : 1 : False : True : Binary
(2, 1) : 0 : None : 1 : False : True : Binary
(2, 2) : 0 : None : 1 : False : True : Binary
(2, 3) : 0 : None : 1 : False : True : Binary
(2, 4) : 0 : None : 1 : False : True : Binary
(2, 5) : 0 : None : 1 : False : True : Binary
(2, 6) : 0 : None : 1 : False : True : Binary
(2, 7) : 0 : None : 1 : False : True : Binary
(2, 8) : 0 : None : 1 : False : True : Binary
(2, 9) : 0 : None : 1 : False : True : Binary
(2, 10) : 0 : None : 1 : False : True : Binary
(3, 0) : 0 : 0.0 : 1 : False : False : Binary
(3, 1) : 0 : 0.0 : 1 : False : False : Binary
(3, 2) : 0 : None : 1 : False : True : Binary
(3, 3) : 0 : None : 1 : False : True : Binary
(3, 4) : 0 : None : 1 : False : True : Binary
(3, 5) : 0 : None : 1 : False : True : Binary
(3, 6) : 0 : None : 1 : False : True : Binary
(3, 7) : 0 : None : 1 : False : True : Binary
(3, 8) : 0 : None : 1 : False : True : Binary
(3, 9) : 0 : None : 1 : False : True : Binary
(3, 10) : 0 : None : 1 : False : True : Binary
(4, 0) : 0 : 0.0 : 1 : False : False : Binary
(4, 1) : 0 : None : 1 : False : True : Binary
(4, 2) : 0 : 0.0 : 1 : False : False : Binary
(4, 3) : 0 : 0.0 : 1 : False : False : Binary
(4, 4) : 0 : None : 1 : False : True : Binary
(4, 5) : 0 : None : 1 : False : True : Binary
(4, 6) : 0 : None : 1 : False : True : Binary
(4, 7) : 0 : None : 1 : False : True : Binary
(4, 8) : 0 : None : 1 : False : True : Binary
(4, 9) : 0 : None : 1 : False : True : Binary
(4, 10) : 0 : None : 1 : False : True : Binary
(5, 0) : 0 : 0.0 : 1 : False : False : Binary
(5, 1) : 0 : 0.0 : 1 : False : False : Binary
(5, 2) : 0 : None : 1 : False : True : Binary
(5, 3) : 0 : 0.0 : 1 : False : False : Binary
(5, 4) : 0 : 0.0 : 1 : False : False : Binary
(5, 5) : 0 : None : 1 : False : True : Binary
(5, 6) : 0 : None : 1 : False : True : Binary
(5, 7) : 0 : None : 1 : False : True : Binary
(5, 8) : 0 : None : 1 : False : True : Binary
(5, 9) : 0 : None : 1 : False : True : Binary
(5, 10) : 0 : None : 1 : False : True : Binary
(6, 0) : 0 : None : 1 : False : True : Binary
(6, 1) : 0 : 0.0 : 1 : False : False : Binary
(6, 2) : 0 : 0.0 : 1 : False : False : Binary
(6, 3) : 0 : 0.0 : 1 : False : False : Binary
(6, 4) : 0 : 0.0 : 1 : False : False : Binary
(6, 5) : 0 : 0.0 : 1 : False : False : Binary
(6, 6) : 0 : None : 1 : False : True : Binary
(6, 7) : 0 : None : 1 : False : True : Binary
(6, 8) : 0 : None : 1 : False : True : Binary
(6, 9) : 0 : None : 1 : False : True : Binary
(6, 10) : 0 : None : 1 : False : True : Binary
(7, 0) : 0 : 0.0 : 1 : False : False : Binary
(7, 1) : 0 : None : 1 : False : True : Binary
(7, 2) : 0 : 0.0 : 1 : False : False : Binary
(7, 3) : 0 : 0.0 : 1 : False : False : Binary
(7, 4) : 0 : 0.0 : 1 : False : False : Binary
(7, 5) : 0 : 0.0 : 1 : False : False : Binary
(7, 6) : 0 : 0.0 : 1 : False : False : Binary
(7, 7) : 0 : None : 1 : False : True : Binary
(7, 8) : 0 : None : 1 : False : True : Binary
(7, 9) : 0 : None : 1 : False : True : Binary
(7, 10) : 0 : None : 1 : False : True : Binary
(8, 0) : 0 : 0.0 : 1 : False : False : Binary
(8, 1) : 0 : 0.0 : 1 : False : False : Binary
(8, 2) : 0 : 0.0 : 1 : False : False : Binary
(8, 3) : 0 : 0.0 : 1 : False : False : Binary
(8, 4) : 0 : 0.0 : 1 : False : False : Binary
(8, 5) : 0 : 0.0 : 1 : False : False : Binary
(8, 6) : 0 : 0.0 : 1 : False : False : Binary
(8, 7) : 0 : 0.0 : 1 : False : False : Binary
(8, 8) : 0 : None : 1 : False : True : Binary
(8, 9) : 0 : None : 1 : False : True : Binary
(8, 10) : 0 : None : 1 : False : True : Binary
(9, 0) : 0 : 0.0 : 1 : False : False : Binary
(9, 1) : 0 : 0.0 : 1 : False : False : Binary
(9, 2) : 0 : 0.0 : 1 : False : False : Binary
(9, 3) : 0 : 0.0 : 1 : False : False : Binary
(9, 4) : 0 : 0.0 : 1 : False : False : Binary
(9, 5) : 0 : 0.0 : 1 : False : False : Binary
(9, 6) : 0 : 0.0 : 1 : False : False : Binary
(9, 7) : 0 : 0.0 : 1 : False : False : Binary
(9, 8) : 0 : 0.0 : 1 : False : False : Binary
(9, 9) : 0 : None : 1 : False : True : Binary
(9, 10) : 0 : None : 1 : False : True : Binary
(10, 0) : 0 : 0.0 : 1 : False : False : Binary
(10, 1) : 0 : 0.0 : 1 : False : False : Binary
(10, 2) : 0 : 0.0 : 1 : False : False : Binary
(10, 3) : 0 : 0.0 : 1 : False : False : Binary
(10, 4) : 0 : 0.0 : 1 : False : False : Binary
(10, 5) : 0 : 0.0 : 1 : False : False : Binary
(10, 6) : 0 : 0.0 : 1 : False : False : Binary
(10, 7) : 0 : 0.0 : 1 : False : False : Binary
(10, 8) : 0 : 0.0 : 1 : False : False : Binary
(10, 9) : 0 : 0.0 : 1 : False : False : Binary
(10, 10) : 0 : None : 1 : False : True : Binary
Objectives:
objective : Size=1, Index=None, Active=True
Key : Active : Value
None : True : 5.0
Constraints:
c1 : Size=47
Key : Lower : Body : Upper
1 : None : -1.0 : 0.0
2 : None : -1.0 : 0.0
3 : None : -1.0 : 0.0
4 : None : -1.0 : 0.0
5 : None : -1.0 : 0.0
6 : None : -1.0 : 0.0
7 : None : -1.0 : 0.0
8 : None : -1.0 : 0.0
9 : None : -1.0 : 0.0
10 : None : -1.0 : 0.0
11 : None : -1.0 : 0.0
12 : None : -1.0 : 0.0
13 : None : -1.0 : 0.0
14 : None : -1.0 : 0.0
15 : None : -1.0 : 0.0
16 : None : -1.0 : 0.0
17 : None : -1.0 : 0.0
18 : None : -1.0 : 0.0
19 : None : -1.0 : 0.0
20 : None : -1.0 : 0.0
21 : None : -1.0 : 0.0
22 : None : -1.0 : 0.0
23 : None : -1.0 : 0.0
24 : None : -1.0 : 0.0
25 : None : -1.0 : 0.0
26 : None : -1.0 : 0.0
27 : None : -1.0 : 0.0
28 : None : -1.0 : 0.0
29 : None : -1.0 : 0.0
30 : None : -1.0 : 0.0
31 : None : -1.0 : 0.0
32 : None : -1.0 : 0.0
33 : None : -1.0 : 0.0
34 : None : -1.0 : 0.0
35 : None : -1.0 : 0.0
36 : None : -1.0 : 0.0
37 : None : -1.0 : 0.0
38 : None : 0.0 : 0.0
39 : None : 0.0 : 0.0
40 : None : 0.0 : 0.0
41 : None : 0.0 : 0.0
42 : None : 0.0 : 0.0
43 : None : 0.0 : 0.0
44 : None : 0.0 : 0.0
45 : None : 0.0 : 0.0
46 : None : 0.0 : 0.0
47 : None : 0.0 : 0.0
c2 : Size=47
Key : Lower : Body : Upper
1 : None : 0.0 : 0.0
2 : None : 0.0 : 0.0
3 : None : 0.0 : 0.0
4 : None : 0.0 : 0.0
5 : None : 0.0 : 0.0
6 : None : 0.0 : 0.0
7 : None : 0.0 : 0.0
8 : None : 0.0 : 0.0
9 : None : 0.0 : 0.0
10 : None : 0.0 : 0.0
11 : None : 0.0 : 0.0
12 : None : 0.0 : 0.0
13 : None : 0.0 : 0.0
14 : None : 0.0 : 0.0
15 : None : 0.0 : 0.0
16 : None : 0.0 : 0.0
17 : None : 0.0 : 0.0
18 : None : 0.0 : 0.0
19 : None : 0.0 : 0.0
20 : None : 0.0 : 0.0
21 : None : 0.0 : 0.0
22 : None : 0.0 : 0.0
23 : None : 0.0 : 0.0
24 : None : 0.0 : 0.0
25 : None : 0.0 : 0.0
26 : None : 0.0 : 0.0
27 : None : 0.0 : 0.0
28 : None : 0.0 : 0.0
29 : None : 0.0 : 0.0
30 : None : 0.0 : 0.0
31 : None : 0.0 : 0.0
32 : None : 0.0 : 0.0
33 : None : 0.0 : 0.0
34 : None : 0.0 : 0.0
35 : None : 0.0 : 0.0
36 : None : 0.0 : 0.0
37 : None : 0.0 : 0.0
38 : None : -1.0 : 0.0
39 : None : -1.0 : 0.0
40 : None : -1.0 : 0.0
41 : None : -1.0 : 0.0
42 : None : -1.0 : 0.0
43 : None : -1.0 : 0.0
44 : None : -1.0 : 0.0
45 : None : -1.0 : 0.0
46 : None : -1.0 : 0.0
47 : None : -1.0 : 0.0
c3 : Size=47
Key : Lower : Body : Upper
1 : None : 0.0 : 0.0
2 : None : 0.0 : 0.0
3 : None : 0.0 : 0.0
4 : None : 0.0 : 0.0
5 : None : 0.0 : 0.0
6 : None : 0.0 : 0.0
7 : None : 0.0 : 0.0
8 : None : 0.0 : 0.0
9 : None : 0.0 : 0.0
10 : None : 0.0 : 0.0
11 : None : 0.0 : 0.0
12 : None : 0.0 : 0.0
13 : None : 0.0 : 0.0
14 : None : 0.0 : 0.0
15 : None : 0.0 : 0.0
16 : None : 0.0 : 0.0
17 : None : 0.0 : 0.0
18 : None : 0.0 : 0.0
19 : None : 0.0 : 0.0
20 : None : 0.0 : 0.0
21 : None : 0.0 : 0.0
22 : None : 0.0 : 0.0
23 : None : 0.0 : 0.0
24 : None : 0.0 : 0.0
25 : None : 0.0 : 0.0
26 : None : 0.0 : 0.0
27 : None : 0.0 : 0.0
28 : None : 0.0 : 0.0
29 : None : 0.0 : 0.0
30 : None : 0.0 : 0.0
31 : None : 0.0 : 0.0
32 : None : 0.0 : 0.0
33 : None : 0.0 : 0.0
34 : None : 0.0 : 0.0
35 : None : 0.0 : 0.0
36 : None : 0.0 : 0.0
37 : None : 0.0 : 0.0
38 : None : 0.0 : 0.0
39 : None : 0.0 : 0.0
40 : None : 0.0 : 0.0
41 : None : 0.0 : 0.0
42 : None : 0.0 : 0.0
43 : None : 0.0 : 0.0
44 : None : 0.0 : 0.0
45 : None : 0.0 : 0.0
46 : None : 0.0 : 0.0
47 : None : 0.0 : 0.0
We observe that the optimal solution of this problem is otherwise, leading to an objective of 5. Notice that this problem has a degenerate optimal solution given that otherwise also leads to the same solution.
Ising model¶
This section introduces the Ising model. We use D-Wave’s dimod and neal packages to define Ising models and solve them with simulated annealing. For integer-programming comparisons, we use Pyomo, an open-source modeling framework with access to linear and nonlinear solvers. The examples use the open-source GLPK solver for mixed-integer linear programming.
Ising problem statement¶
We pose the Ising problem as the following optimization problem:
where we optimize over spins , on a constrained graph , where the quadratic coefficients are and the linear coefficients are . We also include an arbitrary offset of the Ising model .
# These could also be simple lists and numpy matrices
h = {0: 145.0, 1: 122.0, 2: 122.0, 3: 266.0, 4: 266.0, 5: 266.0, 6: 242.5, 7: 266.0, 8: 386.5, 9: 387.0, 10: 386.5}
# Ising coupling matrix (J_{ij})
J = {(0, 3): 24.0, (0, 4): 24.0, (0, 5): 24.0, (0, 7): 24.0, (0, 8): 24.0, (0, 9): 24.0, (0, 10): 24.0, (1, 3): 24.0, (1, 5): 24.0, (1, 6): 24.0, (1, 8): 24.0, (1, 9): 24.0, (1, 10): 24.0, (2, 4): 24.0, (2, 6): 24.0, (2, 7): 24.0, (2, 8): 24.0, (2, 9): 24.0, (2, 10): 24.0, (3, 4): 24.0, (3, 5): 48.0, (3, 6): 24.0, (3, 7): 24.0, (3, 8): 48.0, (3, 9): 48.0, (3, 10): 48.0, (4, 5): 24.0, (4, 6): 24.0, (4, 7): 48.0, (4, 8): 48.0, (4, 9): 48.0, (4, 10): 48.0, (5, 6): 24.0, (5, 7): 24.0, (5, 8): 48.0, (5, 9): 48.0, (5, 10): 48.0, (6, 7): 24.0, (6, 8): 48.0, (6, 9): 48.0, (6, 10): 48.0, (7, 8): 48.0, (7, 9): 48.0, (7, 10): 48.0, (8, 9): 72.0, (8, 10): 72.0, (9, 10): 72.0}
cI = 1319.5
model_ising = dimod.BinaryQuadraticModel.from_ising(h, J, offset=cI)Since the problem is relatively small (11 variables, combinations), we can afford to enumerate all the solutions.
exactSamples = exactSampler.sample(model_ising)plot_enumerate(exactSamples, title='Enumerate all solutions')
plot_energies(exactSamples, title='Enumerate all solutions')
minimum energy: 5.0
minimum energy: 5.0
Let’s now solve this Ising Model via traditional Integer Programming.
# We do not need to worry about the tranformation from Ising to QUBO since dimod takes care of it
Q, c = model_ising.to_qubo()
# Define the model
model_ising_pyo = pyo.ConcreteModel(name='Ising example as an IP, 47-779/785 QuIPML')
I = range(len(h))
J_idx = range(len(h)) # Pyomo index set for Ising couplings
#Define the original variables
model_ising_pyo.x = pyo.Var(I, domain=pyo.Binary)
# Define the edges variables
model_ising_pyo.y = pyo.Var(I, J_idx, domain=pyo.Binary)
obj_expr = c
# add model constraints
model_ising_pyo.c1 = pyo.ConstraintList()
model_ising_pyo.c2 = pyo.ConstraintList()
model_ising_pyo.c3 = pyo.ConstraintList()
for (i,j) in Q.keys():
if i != j:
model_ising_pyo.c1.add(model_ising_pyo.y[i,j] >= model_ising_pyo.x[i] + model_ising_pyo.x[j] - 1)
model_ising_pyo.c2.add(model_ising_pyo.y[i,j] <= model_ising_pyo.x[i])
model_ising_pyo.c3.add(model_ising_pyo.y[i,j] <= model_ising_pyo.x[j])
obj_expr += Q[i,j]*model_ising_pyo.y[i,j]
else:
obj_expr += Q[i,j]*model_ising_pyo.x[i]
# Define the objective function
model_ising_pyo.objective = pyo.Objective(expr = obj_expr, sense=pyo.minimize)
# Print the model
model_ising_pyo.display()Model 'Ising example as an IP, 47-779/785 QuIPML'
Variables:
x : Size=11, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
0 : 0 : None : 1 : False : True : Binary
1 : 0 : None : 1 : False : True : Binary
2 : 0 : None : 1 : False : True : Binary
3 : 0 : None : 1 : False : True : Binary
4 : 0 : None : 1 : False : True : Binary
5 : 0 : None : 1 : False : True : Binary
6 : 0 : None : 1 : False : True : Binary
7 : 0 : None : 1 : False : True : Binary
8 : 0 : None : 1 : False : True : Binary
9 : 0 : None : 1 : False : True : Binary
10 : 0 : None : 1 : False : True : Binary
y : Size=121, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}*{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
(0, 0) : 0 : None : 1 : False : True : Binary
(0, 1) : 0 : None : 1 : False : True : Binary
(0, 2) : 0 : None : 1 : False : True : Binary
(0, 3) : 0 : None : 1 : False : True : Binary
(0, 4) : 0 : None : 1 : False : True : Binary
(0, 5) : 0 : None : 1 : False : True : Binary
(0, 6) : 0 : None : 1 : False : True : Binary
(0, 7) : 0 : None : 1 : False : True : Binary
(0, 8) : 0 : None : 1 : False : True : Binary
(0, 9) : 0 : None : 1 : False : True : Binary
(0, 10) : 0 : None : 1 : False : True : Binary
(1, 0) : 0 : None : 1 : False : True : Binary
(1, 1) : 0 : None : 1 : False : True : Binary
(1, 2) : 0 : None : 1 : False : True : Binary
(1, 3) : 0 : None : 1 : False : True : Binary
(1, 4) : 0 : None : 1 : False : True : Binary
(1, 5) : 0 : None : 1 : False : True : Binary
(1, 6) : 0 : None : 1 : False : True : Binary
(1, 7) : 0 : None : 1 : False : True : Binary
(1, 8) : 0 : None : 1 : False : True : Binary
(1, 9) : 0 : None : 1 : False : True : Binary
(1, 10) : 0 : None : 1 : False : True : Binary
(2, 0) : 0 : None : 1 : False : True : Binary
(2, 1) : 0 : None : 1 : False : True : Binary
(2, 2) : 0 : None : 1 : False : True : Binary
(2, 3) : 0 : None : 1 : False : True : Binary
(2, 4) : 0 : None : 1 : False : True : Binary
(2, 5) : 0 : None : 1 : False : True : Binary
(2, 6) : 0 : None : 1 : False : True : Binary
(2, 7) : 0 : None : 1 : False : True : Binary
(2, 8) : 0 : None : 1 : False : True : Binary
(2, 9) : 0 : None : 1 : False : True : Binary
(2, 10) : 0 : None : 1 : False : True : Binary
(3, 0) : 0 : None : 1 : False : True : Binary
(3, 1) : 0 : None : 1 : False : True : Binary
(3, 2) : 0 : None : 1 : False : True : Binary
(3, 3) : 0 : None : 1 : False : True : Binary
(3, 4) : 0 : None : 1 : False : True : Binary
(3, 5) : 0 : None : 1 : False : True : Binary
(3, 6) : 0 : None : 1 : False : True : Binary
(3, 7) : 0 : None : 1 : False : True : Binary
(3, 8) : 0 : None : 1 : False : True : Binary
(3, 9) : 0 : None : 1 : False : True : Binary
(3, 10) : 0 : None : 1 : False : True : Binary
(4, 0) : 0 : None : 1 : False : True : Binary
(4, 1) : 0 : None : 1 : False : True : Binary
(4, 2) : 0 : None : 1 : False : True : Binary
(4, 3) : 0 : None : 1 : False : True : Binary
(4, 4) : 0 : None : 1 : False : True : Binary
(4, 5) : 0 : None : 1 : False : True : Binary
(4, 6) : 0 : None : 1 : False : True : Binary
(4, 7) : 0 : None : 1 : False : True : Binary
(4, 8) : 0 : None : 1 : False : True : Binary
(4, 9) : 0 : None : 1 : False : True : Binary
(4, 10) : 0 : None : 1 : False : True : Binary
(5, 0) : 0 : None : 1 : False : True : Binary
(5, 1) : 0 : None : 1 : False : True : Binary
(5, 2) : 0 : None : 1 : False : True : Binary
(5, 3) : 0 : None : 1 : False : True : Binary
(5, 4) : 0 : None : 1 : False : True : Binary
(5, 5) : 0 : None : 1 : False : True : Binary
(5, 6) : 0 : None : 1 : False : True : Binary
(5, 7) : 0 : None : 1 : False : True : Binary
(5, 8) : 0 : None : 1 : False : True : Binary
(5, 9) : 0 : None : 1 : False : True : Binary
(5, 10) : 0 : None : 1 : False : True : Binary
(6, 0) : 0 : None : 1 : False : True : Binary
(6, 1) : 0 : None : 1 : False : True : Binary
(6, 2) : 0 : None : 1 : False : True : Binary
(6, 3) : 0 : None : 1 : False : True : Binary
(6, 4) : 0 : None : 1 : False : True : Binary
(6, 5) : 0 : None : 1 : False : True : Binary
(6, 6) : 0 : None : 1 : False : True : Binary
(6, 7) : 0 : None : 1 : False : True : Binary
(6, 8) : 0 : None : 1 : False : True : Binary
(6, 9) : 0 : None : 1 : False : True : Binary
(6, 10) : 0 : None : 1 : False : True : Binary
(7, 0) : 0 : None : 1 : False : True : Binary
(7, 1) : 0 : None : 1 : False : True : Binary
(7, 2) : 0 : None : 1 : False : True : Binary
(7, 3) : 0 : None : 1 : False : True : Binary
(7, 4) : 0 : None : 1 : False : True : Binary
(7, 5) : 0 : None : 1 : False : True : Binary
(7, 6) : 0 : None : 1 : False : True : Binary
(7, 7) : 0 : None : 1 : False : True : Binary
(7, 8) : 0 : None : 1 : False : True : Binary
(7, 9) : 0 : None : 1 : False : True : Binary
(7, 10) : 0 : None : 1 : False : True : Binary
(8, 0) : 0 : None : 1 : False : True : Binary
(8, 1) : 0 : None : 1 : False : True : Binary
(8, 2) : 0 : None : 1 : False : True : Binary
(8, 3) : 0 : None : 1 : False : True : Binary
(8, 4) : 0 : None : 1 : False : True : Binary
(8, 5) : 0 : None : 1 : False : True : Binary
(8, 6) : 0 : None : 1 : False : True : Binary
(8, 7) : 0 : None : 1 : False : True : Binary
(8, 8) : 0 : None : 1 : False : True : Binary
(8, 9) : 0 : None : 1 : False : True : Binary
(8, 10) : 0 : None : 1 : False : True : Binary
(9, 0) : 0 : None : 1 : False : True : Binary
(9, 1) : 0 : None : 1 : False : True : Binary
(9, 2) : 0 : None : 1 : False : True : Binary
(9, 3) : 0 : None : 1 : False : True : Binary
(9, 4) : 0 : None : 1 : False : True : Binary
(9, 5) : 0 : None : 1 : False : True : Binary
(9, 6) : 0 : None : 1 : False : True : Binary
(9, 7) : 0 : None : 1 : False : True : Binary
(9, 8) : 0 : None : 1 : False : True : Binary
(9, 9) : 0 : None : 1 : False : True : Binary
(9, 10) : 0 : None : 1 : False : True : Binary
(10, 0) : 0 : None : 1 : False : True : Binary
(10, 1) : 0 : None : 1 : False : True : Binary
(10, 2) : 0 : None : 1 : False : True : Binary
(10, 3) : 0 : None : 1 : False : True : Binary
(10, 4) : 0 : None : 1 : False : True : Binary
(10, 5) : 0 : None : 1 : False : True : Binary
(10, 6) : 0 : None : 1 : False : True : Binary
(10, 7) : 0 : None : 1 : False : True : Binary
(10, 8) : 0 : None : 1 : False : True : Binary
(10, 9) : 0 : None : 1 : False : True : Binary
(10, 10) : 0 : None : 1 : False : True : Binary
Objectives:
objective : Size=1, Index=None, Active=True
ERROR: evaluating object as numeric value: y[3,0]
(object: <class 'pyomo.core.base.var.VarData'>)
No value for uninitialized VarData object y[3,0]
ERROR: evaluating object as numeric value: objective
(object: <class 'pyomo.core.base.objective.ScalarObjective'>)
No value for uninitialized VarData object y[3,0]
Key : Active : Value
None : None : None
Constraints:
c1 : Size=47
Key : Lower : Body : Upper
1 : None : None : 0.0
2 : None : None : 0.0
3 : None : None : 0.0
4 : None : None : 0.0
5 : None : None : 0.0
6 : None : None : 0.0
7 : None : None : 0.0
8 : None : None : 0.0
9 : None : None : 0.0
10 : None : None : 0.0
11 : None : None : 0.0
12 : None : None : 0.0
13 : None : None : 0.0
14 : None : None : 0.0
15 : None : None : 0.0
16 : None : None : 0.0
17 : None : None : 0.0
18 : None : None : 0.0
19 : None : None : 0.0
20 : None : None : 0.0
21 : None : None : 0.0
22 : None : None : 0.0
23 : None : None : 0.0
24 : None : None : 0.0
25 : None : None : 0.0
26 : None : None : 0.0
27 : None : None : 0.0
28 : None : None : 0.0
29 : None : None : 0.0
30 : None : None : 0.0
31 : None : None : 0.0
32 : None : None : 0.0
33 : None : None : 0.0
34 : None : None : 0.0
35 : None : None : 0.0
36 : None : None : 0.0
37 : None : None : 0.0
38 : None : None : 0.0
39 : None : None : 0.0
40 : None : None : 0.0
41 : None : None : 0.0
42 : None : None : 0.0
43 : None : None : 0.0
44 : None : None : 0.0
45 : None : None : 0.0
46 : None : None : 0.0
47 : None : None : 0.0
c2 : Size=47
Key : Lower : Body : Upper
1 : None : None : 0.0
2 : None : None : 0.0
3 : None : None : 0.0
4 : None : None : 0.0
5 : None : None : 0.0
6 : None : None : 0.0
7 : None : None : 0.0
8 : None : None : 0.0
9 : None : None : 0.0
10 : None : None : 0.0
11 : None : None : 0.0
12 : None : None : 0.0
13 : None : None : 0.0
14 : None : None : 0.0
15 : None : None : 0.0
16 : None : None : 0.0
17 : None : None : 0.0
18 : None : None : 0.0
19 : None : None : 0.0
20 : None : None : 0.0
21 : None : None : 0.0
22 : None : None : 0.0
23 : None : None : 0.0
24 : None : None : 0.0
25 : None : None : 0.0
26 : None : None : 0.0
27 : None : None : 0.0
28 : None : None : 0.0
29 : None : None : 0.0
30 : None : None : 0.0
31 : None : None : 0.0
32 : None : None : 0.0
33 : None : None : 0.0
34 : None : None : 0.0
35 : None : None : 0.0
36 : None : None : 0.0
37 : None : None : 0.0
38 : None : None : 0.0
39 : None : None : 0.0
40 : None : None : 0.0
41 : None : None : 0.0
42 : None : None : 0.0
43 : None : None : 0.0
44 : None : None : 0.0
45 : None : None : 0.0
46 : None : None : 0.0
47 : None : None : 0.0
c3 : Size=47
Key : Lower : Body : Upper
1 : None : None : 0.0
2 : None : None : 0.0
3 : None : None : 0.0
4 : None : None : 0.0
5 : None : None : 0.0
6 : None : None : 0.0
7 : None : None : 0.0
8 : None : None : 0.0
9 : None : None : 0.0
10 : None : None : 0.0
11 : None : None : 0.0
12 : None : None : 0.0
13 : None : None : 0.0
14 : None : None : 0.0
15 : None : None : 0.0
16 : None : None : 0.0
17 : None : None : 0.0
18 : None : None : 0.0
19 : None : None : 0.0
20 : None : None : 0.0
21 : None : None : 0.0
22 : None : None : 0.0
23 : None : None : 0.0
24 : None : None : 0.0
25 : None : None : 0.0
26 : None : None : 0.0
27 : None : None : 0.0
28 : None : None : 0.0
29 : None : None : 0.0
30 : None : None : 0.0
31 : None : None : 0.0
32 : None : None : 0.0
33 : None : None : 0.0
34 : None : None : 0.0
35 : None : None : 0.0
36 : None : None : 0.0
37 : None : None : 0.0
38 : None : None : 0.0
39 : None : None : 0.0
40 : None : None : 0.0
41 : None : None : 0.0
42 : None : None : 0.0
43 : None : None : 0.0
44 : None : None : 0.0
45 : None : None : 0.0
46 : None : None : 0.0
47 : None : None : 0.0
# We obtain the solution with GLPK
result_obj = opt_glpk.solve(model_ising_pyo, tee=False)
model_ising_pyo.display()Model 'Ising example as an IP, 47-779/785 QuIPML'
Variables:
x : Size=11, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
0 : 0 : 0.0 : 1 : False : False : Binary
1 : 0 : 0.0 : 1 : False : False : Binary
2 : 0 : 0.0 : 1 : False : False : Binary
3 : 0 : 0.0 : 1 : False : False : Binary
4 : 0 : 0.0 : 1 : False : False : Binary
5 : 0 : 0.0 : 1 : False : False : Binary
6 : 0 : 0.0 : 1 : False : False : Binary
7 : 0 : 0.0 : 1 : False : False : Binary
8 : 0 : 0.0 : 1 : False : False : Binary
9 : 0 : 0.0 : 1 : False : False : Binary
10 : 0 : 1.0 : 1 : False : False : Binary
y : Size=121, Index={0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}*{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
Key : Lower : Value : Upper : Fixed : Stale : Domain
(0, 0) : 0 : None : 1 : False : True : Binary
(0, 1) : 0 : None : 1 : False : True : Binary
(0, 2) : 0 : None : 1 : False : True : Binary
(0, 3) : 0 : None : 1 : False : True : Binary
(0, 4) : 0 : None : 1 : False : True : Binary
(0, 5) : 0 : None : 1 : False : True : Binary
(0, 6) : 0 : None : 1 : False : True : Binary
(0, 7) : 0 : None : 1 : False : True : Binary
(0, 8) : 0 : None : 1 : False : True : Binary
(0, 9) : 0 : None : 1 : False : True : Binary
(0, 10) : 0 : None : 1 : False : True : Binary
(1, 0) : 0 : None : 1 : False : True : Binary
(1, 1) : 0 : None : 1 : False : True : Binary
(1, 2) : 0 : None : 1 : False : True : Binary
(1, 3) : 0 : 0.0 : 1 : False : False : Binary
(1, 4) : 0 : None : 1 : False : True : Binary
(1, 5) : 0 : 0.0 : 1 : False : False : Binary
(1, 6) : 0 : None : 1 : False : True : Binary
(1, 7) : 0 : None : 1 : False : True : Binary
(1, 8) : 0 : 0.0 : 1 : False : False : Binary
(1, 9) : 0 : 0.0 : 1 : False : False : Binary
(1, 10) : 0 : 0.0 : 1 : False : False : Binary
(2, 0) : 0 : None : 1 : False : True : Binary
(2, 1) : 0 : None : 1 : False : True : Binary
(2, 2) : 0 : None : 1 : False : True : Binary
(2, 3) : 0 : None : 1 : False : True : Binary
(2, 4) : 0 : 0.0 : 1 : False : False : Binary
(2, 5) : 0 : None : 1 : False : True : Binary
(2, 6) : 0 : 0.0 : 1 : False : False : Binary
(2, 7) : 0 : 0.0 : 1 : False : False : Binary
(2, 8) : 0 : 0.0 : 1 : False : False : Binary
(2, 9) : 0 : 0.0 : 1 : False : False : Binary
(2, 10) : 0 : 0.0 : 1 : False : False : Binary
(3, 0) : 0 : 0.0 : 1 : False : False : Binary
(3, 1) : 0 : None : 1 : False : True : Binary
(3, 2) : 0 : None : 1 : False : True : Binary
(3, 3) : 0 : None : 1 : False : True : Binary
(3, 4) : 0 : None : 1 : False : True : Binary
(3, 5) : 0 : None : 1 : False : True : Binary
(3, 6) : 0 : None : 1 : False : True : Binary
(3, 7) : 0 : None : 1 : False : True : Binary
(3, 8) : 0 : None : 1 : False : True : Binary
(3, 9) : 0 : None : 1 : False : True : Binary
(3, 10) : 0 : None : 1 : False : True : Binary
(4, 0) : 0 : 0.0 : 1 : False : False : Binary
(4, 1) : 0 : None : 1 : False : True : Binary
(4, 2) : 0 : None : 1 : False : True : Binary
(4, 3) : 0 : 0.0 : 1 : False : False : Binary
(4, 4) : 0 : None : 1 : False : True : Binary
(4, 5) : 0 : None : 1 : False : True : Binary
(4, 6) : 0 : None : 1 : False : True : Binary
(4, 7) : 0 : None : 1 : False : True : Binary
(4, 8) : 0 : None : 1 : False : True : Binary
(4, 9) : 0 : None : 1 : False : True : Binary
(4, 10) : 0 : None : 1 : False : True : Binary
(5, 0) : 0 : 0.0 : 1 : False : False : Binary
(5, 1) : 0 : None : 1 : False : True : Binary
(5, 2) : 0 : None : 1 : False : True : Binary
(5, 3) : 0 : 0.0 : 1 : False : False : Binary
(5, 4) : 0 : 0.0 : 1 : False : False : Binary
(5, 5) : 0 : None : 1 : False : True : Binary
(5, 6) : 0 : None : 1 : False : True : Binary
(5, 7) : 0 : None : 1 : False : True : Binary
(5, 8) : 0 : None : 1 : False : True : Binary
(5, 9) : 0 : None : 1 : False : True : Binary
(5, 10) : 0 : None : 1 : False : True : Binary
(6, 0) : 0 : None : 1 : False : True : Binary
(6, 1) : 0 : 0.0 : 1 : False : False : Binary
(6, 2) : 0 : None : 1 : False : True : Binary
(6, 3) : 0 : 0.0 : 1 : False : False : Binary
(6, 4) : 0 : 0.0 : 1 : False : False : Binary
(6, 5) : 0 : 0.0 : 1 : False : False : Binary
(6, 6) : 0 : None : 1 : False : True : Binary
(6, 7) : 0 : 0.0 : 1 : False : False : Binary
(6, 8) : 0 : 0.0 : 1 : False : False : Binary
(6, 9) : 0 : 0.0 : 1 : False : False : Binary
(6, 10) : 0 : 0.0 : 1 : False : False : Binary
(7, 0) : 0 : 0.0 : 1 : False : False : Binary
(7, 1) : 0 : None : 1 : False : True : Binary
(7, 2) : 0 : None : 1 : False : True : Binary
(7, 3) : 0 : 0.0 : 1 : False : False : Binary
(7, 4) : 0 : 0.0 : 1 : False : False : Binary
(7, 5) : 0 : 0.0 : 1 : False : False : Binary
(7, 6) : 0 : None : 1 : False : True : Binary
(7, 7) : 0 : None : 1 : False : True : Binary
(7, 8) : 0 : None : 1 : False : True : Binary
(7, 9) : 0 : None : 1 : False : True : Binary
(7, 10) : 0 : None : 1 : False : True : Binary
(8, 0) : 0 : 0.0 : 1 : False : False : Binary
(8, 1) : 0 : None : 1 : False : True : Binary
(8, 2) : 0 : None : 1 : False : True : Binary
(8, 3) : 0 : 0.0 : 1 : False : False : Binary
(8, 4) : 0 : 0.0 : 1 : False : False : Binary
(8, 5) : 0 : 0.0 : 1 : False : False : Binary
(8, 6) : 0 : None : 1 : False : True : Binary
(8, 7) : 0 : 0.0 : 1 : False : False : Binary
(8, 8) : 0 : None : 1 : False : True : Binary
(8, 9) : 0 : None : 1 : False : True : Binary
(8, 10) : 0 : None : 1 : False : True : Binary
(9, 0) : 0 : 0.0 : 1 : False : False : Binary
(9, 1) : 0 : None : 1 : False : True : Binary
(9, 2) : 0 : None : 1 : False : True : Binary
(9, 3) : 0 : 0.0 : 1 : False : False : Binary
(9, 4) : 0 : 0.0 : 1 : False : False : Binary
(9, 5) : 0 : 0.0 : 1 : False : False : Binary
(9, 6) : 0 : None : 1 : False : True : Binary
(9, 7) : 0 : 0.0 : 1 : False : False : Binary
(9, 8) : 0 : 0.0 : 1 : False : False : Binary
(9, 9) : 0 : None : 1 : False : True : Binary
(9, 10) : 0 : None : 1 : False : True : Binary
(10, 0) : 0 : 0.0 : 1 : False : False : Binary
(10, 1) : 0 : None : 1 : False : True : Binary
(10, 2) : 0 : None : 1 : False : True : Binary
(10, 3) : 0 : 0.0 : 1 : False : False : Binary
(10, 4) : 0 : 0.0 : 1 : False : False : Binary
(10, 5) : 0 : 0.0 : 1 : False : False : Binary
(10, 6) : 0 : None : 1 : False : True : Binary
(10, 7) : 0 : 0.0 : 1 : False : False : Binary
(10, 8) : 0 : 0.0 : 1 : False : False : Binary
(10, 9) : 0 : 0.0 : 1 : False : False : Binary
(10, 10) : 0 : None : 1 : False : True : Binary
Objectives:
objective : Size=1, Index=None, Active=True
Key : Active : Value
None : True : 5.0
Constraints:
c1 : Size=47
Key : Lower : Body : Upper
1 : None : -1.0 : 0.0
2 : None : -1.0 : 0.0
3 : None : -1.0 : 0.0
4 : None : -1.0 : 0.0
5 : None : -1.0 : 0.0
6 : None : -1.0 : 0.0
7 : None : -1.0 : 0.0
8 : None : -1.0 : 0.0
9 : None : -1.0 : 0.0
10 : None : -1.0 : 0.0
11 : None : -1.0 : 0.0
12 : None : -1.0 : 0.0
13 : None : -1.0 : 0.0
14 : None : -1.0 : 0.0
15 : None : -1.0 : 0.0
16 : None : -1.0 : 0.0
17 : None : -1.0 : 0.0
18 : None : -1.0 : 0.0
19 : None : -1.0 : 0.0
20 : None : -1.0 : 0.0
21 : None : -1.0 : 0.0
22 : None : 0.0 : 0.0
23 : None : 0.0 : 0.0
24 : None : 0.0 : 0.0
25 : None : 0.0 : 0.0
26 : None : 0.0 : 0.0
27 : None : 0.0 : 0.0
28 : None : 0.0 : 0.0
29 : None : -1.0 : 0.0
30 : None : -1.0 : 0.0
31 : None : -1.0 : 0.0
32 : None : -1.0 : 0.0
33 : None : 0.0 : 0.0
34 : None : -1.0 : 0.0
35 : None : -1.0 : 0.0
36 : None : -1.0 : 0.0
37 : None : -1.0 : 0.0
38 : None : -1.0 : 0.0
39 : None : -1.0 : 0.0
40 : None : 0.0 : 0.0
41 : None : -1.0 : 0.0
42 : None : -1.0 : 0.0
43 : None : -1.0 : 0.0
44 : None : -1.0 : 0.0
45 : None : -1.0 : 0.0
46 : None : 0.0 : 0.0
47 : None : -1.0 : 0.0
c2 : Size=47
Key : Lower : Body : Upper
1 : None : 0.0 : 0.0
2 : None : 0.0 : 0.0
3 : None : 0.0 : 0.0
4 : None : 0.0 : 0.0
5 : None : 0.0 : 0.0
6 : None : 0.0 : 0.0
7 : None : 0.0 : 0.0
8 : None : 0.0 : 0.0
9 : None : 0.0 : 0.0
10 : None : 0.0 : 0.0
11 : None : 0.0 : 0.0
12 : None : 0.0 : 0.0
13 : None : 0.0 : 0.0
14 : None : 0.0 : 0.0
15 : None : 0.0 : 0.0
16 : None : 0.0 : 0.0
17 : None : 0.0 : 0.0
18 : None : 0.0 : 0.0
19 : None : 0.0 : 0.0
20 : None : 0.0 : 0.0
21 : None : 0.0 : 0.0
22 : None : -1.0 : 0.0
23 : None : -1.0 : 0.0
24 : None : -1.0 : 0.0
25 : None : -1.0 : 0.0
26 : None : -1.0 : 0.0
27 : None : -1.0 : 0.0
28 : None : -1.0 : 0.0
29 : None : 0.0 : 0.0
30 : None : 0.0 : 0.0
31 : None : 0.0 : 0.0
32 : None : 0.0 : 0.0
33 : None : 0.0 : 0.0
34 : None : 0.0 : 0.0
35 : None : 0.0 : 0.0
36 : None : 0.0 : 0.0
37 : None : 0.0 : 0.0
38 : None : 0.0 : 0.0
39 : None : 0.0 : 0.0
40 : None : 0.0 : 0.0
41 : None : 0.0 : 0.0
42 : None : 0.0 : 0.0
43 : None : 0.0 : 0.0
44 : None : 0.0 : 0.0
45 : None : 0.0 : 0.0
46 : None : 0.0 : 0.0
47 : None : 0.0 : 0.0
c3 : Size=47
Key : Lower : Body : Upper
1 : None : 0.0 : 0.0
2 : None : 0.0 : 0.0
3 : None : 0.0 : 0.0
4 : None : 0.0 : 0.0
5 : None : 0.0 : 0.0
6 : None : 0.0 : 0.0
7 : None : 0.0 : 0.0
8 : None : 0.0 : 0.0
9 : None : 0.0 : 0.0
10 : None : 0.0 : 0.0
11 : None : 0.0 : 0.0
12 : None : 0.0 : 0.0
13 : None : 0.0 : 0.0
14 : None : 0.0 : 0.0
15 : None : 0.0 : 0.0
16 : None : 0.0 : 0.0
17 : None : 0.0 : 0.0
18 : None : 0.0 : 0.0
19 : None : 0.0 : 0.0
20 : None : 0.0 : 0.0
21 : None : 0.0 : 0.0
22 : None : 0.0 : 0.0
23 : None : 0.0 : 0.0
24 : None : 0.0 : 0.0
25 : None : 0.0 : 0.0
26 : None : 0.0 : 0.0
27 : None : 0.0 : 0.0
28 : None : 0.0 : 0.0
29 : None : 0.0 : 0.0
30 : None : 0.0 : 0.0
31 : None : 0.0 : 0.0
32 : None : 0.0 : 0.0
33 : None : -1.0 : 0.0
34 : None : 0.0 : 0.0
35 : None : 0.0 : 0.0
36 : None : 0.0 : 0.0
37 : None : 0.0 : 0.0
38 : None : 0.0 : 0.0
39 : None : 0.0 : 0.0
40 : None : -1.0 : 0.0
41 : None : 0.0 : 0.0
42 : None : 0.0 : 0.0
43 : None : 0.0 : 0.0
44 : None : 0.0 : 0.0
45 : None : 0.0 : 0.0
46 : None : -1.0 : 0.0
47 : None : 0.0 : 0.0
We observe that the optimal solution of this problem is otherwise, leading to an objective of 5. Notice that this problem has a degenerate optimal solution given that otherwise also leads to the same solution.
From QUBO to Ising¶
QUBO and Ising models describe the same binary search space with different variable conventions. The QUBO form uses variables in {0, 1}, while the Ising form uses spins in {-1, +1}. Moving between them changes the linear terms, quadratic terms, and offset, but not the underlying set of assignments. This equivalence lets us compare exact integer-programming solves with simulated annealing on the same problem.
We can also solve this problem using Simulated Annealing
simAnnSampler = neal.SimulatedAnnealingSampler()
# A fixed seed makes the local sampling repeatable.
simAnnSamples = simAnnSampler.sample(model_ising, num_reads=1000, seed=314159)plot_enumerate(simAnnSamples, title='Simulated annealing in default parameters')
plot_energies(simAnnSamples, title='Simulated annealing in default parameters')
minimum energy: 5.0
minimum energy: 5.0
simAnnSamples.info{'beta_range': [np.float64(0.00041111932417553103),
np.float64(0.1458971970580513)],
'beta_schedule_type': 'geometric',
'timing': {'preprocessing_ns': 1258757,
'sampling_ns': 174347111,
'postprocessing_ns': 293102}}From penalty models to graph coloring¶
The previous examples used penalties to enforce algebraic constraints. Graph coloring uses the same idea with one-hot variables: each vertex must choose exactly one color, and adjacent vertices may not share a color. Violating either rule adds a penalty to the QUBO objective, so low-energy samples correspond to valid colorings when the penalties are large enough.
Let’s solve the graph coloring problem using QUBO.
Vertex -coloring of graphs¶
Given a graph , where is the set of vertices and is the set of edges of , and a positive integer , we ask if it is possible to assign a color to every vertex from , such that adjacent vertices have different colors assigned.
has 12 vertices and 23 edges. We ask if the graph is 3–colorable. Let’s first encode and using Python data structures:
Note: This tutorial is heavily inspired in D-Wave’s Map coloring of Canada found here.
# The graph coloring QUBO below uses dimod directly; no CSP package is required.V = range(1, 12+1)
E = [(1,2),(2,3),(1,4),(1,6),(1,12),(2,5),(2,7),(3,8),(3,10),(4,11),(4,9),(5,6),(6,7),(7,8),(8,9),(9,10),(10,11),(11,12),(5,12),(5,9),(6,10),(7,11),(8,12)]
layout = {i: [np.cos((2*i+1)*np.pi/8),np.sin((2*i+1)*np.pi/8)] for i in np.arange(5,13)}
layout[1] = [-1.5,1.5]
layout[2] = [1.5,1.5]
layout[3] = [1.5,-1.5]
layout[4] = [-1.5,-1.5]
G = nx.Graph()
G.add_edges_from(E)
nx.draw(G, with_labels=True, pos=layout)
colors = 3
exactly_one_penalty = 1.0
edge_penalty = 1.0
def color_var(node, color):
return f'x{node},{color}'
def build_graph_coloring_bqm(
nodes,
edges,
colors,
exactly_one_penalty=1.0,
edge_penalty=1.0,
):
bqm = dimod.BinaryQuadraticModel({}, {}, 0.0, dimod.BINARY)
# Add A * (1 - sum_c x_{v,c})^2 for each node v.
for node in nodes:
variables = [color_var(node, color) for color in range(colors)]
bqm.offset += exactly_one_penalty
for variable in variables:
bqm.add_linear(variable, -exactly_one_penalty)
for first in range(colors):
for second in range(first + 1, colors):
bqm.add_quadratic(
variables[first],
variables[second],
2 * exactly_one_penalty,
)
# Add B * x_{u,c} * x_{v,c} for each edge (u, v) and color c.
for v, u in edges:
for color in range(colors):
bqm.add_quadratic(color_var(v, color), color_var(u, color), edge_penalty)
return bqm
def is_valid_coloring(sample):
for node in V:
selected = [color for color in range(colors) if sample[color_var(node, color)]]
if len(selected) != 1:
return False
for v, u in E:
for color in range(colors):
if sample[color_var(v, color)] and sample[color_var(u, color)]:
return False
return TrueDefining the Binary Quadratic Model (QUBO) directly from the graph coloring penalties we have:
bqm = build_graph_coloring_bqm(V, E, colors, exactly_one_penalty, edge_penalty)
simAnnSamples = simAnnSampler.sample(bqm, num_reads=2000, num_sweeps=1000, seed=314159)plot_enumerate(simAnnSamples, title='Simulated annealing in default parameters')
plot_energies(simAnnSamples, title='Simulated annealing in default parameters')
minimum energy: 0.0
minimum energy: 0.0
With this direct penalty model, every feasible 3-coloring has energy 0. Any sample with positive energy violates at least one exactly-one-color or adjacent-color constraint.
# Check that a good solution was found
sample = None
for datum in simAnnSamples.data(['sample'], sorted_by='energy'):
if is_valid_coloring(datum.sample):
sample = datum.sample
break
if sample is None:
raise RuntimeError("Failed to color map. Try sampling again.")
else:
print(sample){'x1,0': np.int8(0), 'x1,1': np.int8(1), 'x1,2': np.int8(0), 'x10,0': np.int8(0), 'x10,1': np.int8(1), 'x10,2': np.int8(0), 'x11,0': np.int8(1), 'x11,1': np.int8(0), 'x11,2': np.int8(0), 'x12,0': np.int8(0), 'x12,1': np.int8(0), 'x12,2': np.int8(1), 'x2,0': np.int8(1), 'x2,1': np.int8(0), 'x2,2': np.int8(0), 'x3,0': np.int8(0), 'x3,1': np.int8(0), 'x3,2': np.int8(1), 'x4,0': np.int8(0), 'x4,1': np.int8(0), 'x4,2': np.int8(1), 'x5,0': np.int8(0), 'x5,1': np.int8(1), 'x5,2': np.int8(0), 'x6,0': np.int8(1), 'x6,1': np.int8(0), 'x6,2': np.int8(0), 'x7,0': np.int8(0), 'x7,1': np.int8(0), 'x7,2': np.int8(1), 'x8,0': np.int8(0), 'x8,1': np.int8(1), 'x8,2': np.int8(0), 'x9,0': np.int8(1), 'x9,1': np.int8(0), 'x9,2': np.int8(0)}
# Function that plots a returned sample
def plot_map(sample):
# Translate from binary to integer color representation
color_map = {}
for node in V:
for i in range(colors):
if sample[color_var(node, i)]:
color_map[node] = i
# Plot the sample with color-coded nodes
node_colors = [color_map.get(node) for node in G.nodes()]
nx.draw(G, with_labels=True, pos=layout, node_color=node_colors)
plt.show()plot_map(sample)
Practice checkpoints¶
Use these checkpoints during the workshop to test the main ideas before moving on.
# =============================================================================
# EXERCISE 1: Build a small QUBO by hand
# =============================================================================
# Write the QUBO coefficients for a three-variable penalty model before using any helper function. Check that the lowest-energy assignment satisfies the original constraint.
#
# Hint: Enumerate all eight assignments to validate the coefficients.
#
# Your code here:
# =============================================================================Notebook Cell
# SOLUTION (hidden in workshop version):
# Build and enumerate a three-variable one-hot QUBO by hand.
import itertools
linear = {0: -1.0, 1: -1.0, 2: -1.0}
quadratic = {(0, 1): 2.0, (0, 2): 2.0, (1, 2): 2.0}
offset = 1.0
def energy(bits):
return (
offset
+ sum(linear[i] * bits[i] for i in linear)
+ sum(weight * bits[i] * bits[j] for (i, j), weight in quadratic.items())
)
best = min(
((bits, energy(bits)) for bits in itertools.product([0, 1], repeat=3)),
key=lambda item: (item[1], -item[0][0], -item[0][1], -item[0][2]),
)
feasible = sum(best[0]) == 1
print(f"best_assignment = {best[0]}; energy = {best[1]}; feasible = {str(feasible).lower()}")
best_assignment = (1, 0, 0); energy = 0.0; feasible = true
# =============================================================================
# EXERCISE 2: Stress-test the penalty parameter
# =============================================================================
# Lower rho below the sufficient bound and sample or enumerate the resulting BQM. Record the first infeasible assignment that becomes competitive with a feasible one.
#
# Hint: Keep the objective and constraints fixed so rho is the only changed input.
#
# Your code here:
# =============================================================================Notebook Cell
# SOLUTION (hidden in workshop version):
# Lower rho on a tiny equality-constrained binary model and detect infeasibility.
import itertools
import numpy as np
A_small = np.array([[1, 1, 0], [0, 1, 1]])
b_small = np.array([1, 1])
c_small = np.array([1.0, 2.0, 1.0])
def penalized_energy(bits, rho_value):
x = np.array(bits)
residual = A_small @ x - b_small
return float(c_small @ x + rho_value * (residual @ residual))
rho_value = 0.25
best = min(
(
(bits, penalized_energy(bits, rho_value), np.array_equal(A_small @ np.array(bits), b_small))
for bits in itertools.product([0, 1], repeat=3)
),
key=lambda item: item[1],
)
print(f"rho = {rho_value}; best_assignment = {best[0]}; energy = {best[1]}; feasible = {str(best[2]).lower()}")
rho = 0.25; best_assignment = (0, 0, 0); energy = 0.5; feasible = false
# =============================================================================
# EXERCISE 3: Extend graph coloring
# =============================================================================
# Change the graph-coloring example to use three colors. Update the variable naming and validity checks so exactly one color is assigned per node.
#
# Hint: The exactly-one penalties scale with the number of colors, not the number of nodes.
#
# Your code here:
# =============================================================================Notebook Cell
# SOLUTION (hidden in workshop version):
# Create variables and constraints for a three-color graph-coloring BQM sketch.
nodes = [1, 2, 3]
edges = [(1, 2), (2, 3)]
colors = ["red", "green", "blue"]
variables = {(node, color): f"x_{node}_{color}" for node in nodes for color in colors}
one_hot_groups = [[variables[node, color] for color in colors] for node in nodes]
edge_conflicts = [
(variables[u, color], variables[v, color])
for u, v in edges
for color in colors
]
format_group = lambda group: ",".join(group)
format_conflict = lambda conflict: f"{conflict[0]}-{conflict[1]}"
print(
"one_hot_groups = "
+ " | ".join(format_group(group) for group in one_hot_groups)
+ "; edge_conflicts = "
+ " | ".join(format_conflict(conflict) for conflict in edge_conflicts)
)
one_hot_groups = x_1_red,x_1_green,x_1_blue | x_2_red,x_2_green,x_2_blue | x_3_red,x_3_green,x_3_blue; edge_conflicts = x_1_red-x_2_red | x_1_green-x_2_green | x_1_blue-x_2_blue | x_2_red-x_3_red | x_2_green-x_3_green | x_2_blue-x_3_blue
The next lecture develops classical augmentation methods. We return to this QUBO for quantum annealing and an optional QAOA comparison in Lecture 4: D-Wave and quantum methods.
Summary¶
In this notebook we:
Defined QUBO and Ising objectives and identified how coefficients encode binary optimization problems.
Reformulated constrained binary problems with penalty terms so they can be sampled as unconstrained BQMs.
Used
dimodand simulated annealing to sample candidate solutions and inspect their energies.Modeled graph coloring directly as a BQM and validated feasible color assignments without a CSP helper package.
Learning objectives met: You practiced reading QUBO/Ising forms, building BQMs, applying penalty reformulations, and checking sampled solutions against the original constraints.
Next steps: Proceed to Notebook 3: GAMA to compare QUBO-style sampling with Graver-basis augmentation for structured integer programs.
Further reading:
David E. Bernal Neira — Davidson School of Chemical Engineering, Purdue University
Pedro Maciel Xavier — Davidson School of Chemical Engineering, Purdue University
Benjamin J. L. Murray — Davidson School of Chemical Engineering, Purdue University; Undergraduate Research Assistant
Acknowledgments¶
This notebook was developed by: