Formulation BenchmarksΒΆ
We have benchmarked the formulation of mathematical optimization models provided as a library in Python against the Amplify SDK. We measured the execution time to create a model and output it as QUBO using the traveling salesperson problem formulation as an example, assuming using the QUBO solver.
However, each library, including the Amplify SDK, covers different features depending on the formulation method. Here is a summary of the features and policies of each library, formulated as follows.
- |
dimod |
dimod |
dimod |
||
---|---|---|---|---|---|
Symbolic operation |
β |
β |
β |
β |
β |
Objective function |
β |
β |
β |
β |
β |
Constraint |
β |
β ** |
β*** |
β*** |
β * |
Higher order polynomial |
β |
β |
β |
β |
β |
Coefficient matrix |
β |
β |
β |
β |
β |
Variable type |
B/S/I/R |
B/S |
B/S |
B/S |
B/S/I/R |
Supported machines |
Various |
Depends on user |
D-Wave only |
D-Wave only |
D-Wave only |
B: Binary, S: Ising Spin, I: Integer, R: Real
*: Model creation only, as QUBO output (conversion of constraints to the penalty functions) is not available.
**: Must define penalty function
***: Constraints are expressed by adding penalty functions to the objective function
Amplify
import amplify
def tsp_for_amplify(ncity: int, distances: np.ndarray, dmax: float):
q = amplify.VariableGenerator().array("Binary", ncity + 1, ncity)
q[-1, :] = q[0, :]
# Objective function
objective: amplify.Poly = amplify.einsum(
"ij,ki,kj->", distances, q[:-1], q[1:]
)
# Constraints
constraints: amplify.ConstraintList = amplify.one_hot(
q[:-1], axis=1
) + amplify.one_hot(q[:-1], axis=0)
return objective + dmax * constraints
class BenchTspAmplify:
def create_model(self, ncity: int, distances: np.ndarray, dmax: float):
self.model = tsp_for_amplify(ncity, distances, dmax)
def to_qubo(self):
self.model.to_unconstrained_poly()
PyQUBO
import pyqubo
def tsp_for_pyqubo(ncity: int, distances: np.ndarray, dmax: float):
# from https://github.com/recruit-communications/pyqubo/blob/master/notebooks/TSP.ipynb
# NOTE: https://github.com/recruit-communications/pyqubo/blob/master/benchmark/benchmark.py
# is not valid for TSP
x = pyqubo.Array.create("c", (ncity, ncity), "BINARY")
# Constraint not to visit more than two cities at the same time.
time_const = 0.0
for i in range(ncity):
# If you wrap the hamiltonian by Const(...), this part is recognized as constraint
time_const += pyqubo.Constraint(
(sum(x[i, j] for j in range(ncity)) - 1) ** 2, label=f"time{i}"
)
# Constraint not to visit the same city more than twice.
city_const = 0.0
for j in range(ncity):
city_const += pyqubo.Constraint(
(sum(x[i, j] for i in range(ncity)) - 1) ** 2, label=f"city{j}"
)
# distance of route
feed_dict = {}
distance = 0.0
for i in range(ncity):
for j in range(ncity):
for k in range(ncity):
# we set the constant distance
distance += distances[i, j] * x[k, i] * x[(k + 1) % ncity, j]
# Construct hamiltonian
A = pyqubo.Placeholder("A")
H = distance + A * (time_const + city_const)
feed_dict["A"] = dmax
# Compile model
return H.compile(), feed_dict
class BenchTspPyQubo:
def create_model(self, ncity: int, distances: np.ndarray, dmax: float):
self.model, self._feed_dict = tsp_for_pyqubo(ncity, distances, dmax)
def to_qubo(self):
self.model.to_qubo(index_label=False, feed_dict=self._feed_dict)
dimod BQM (index)
import dimod
def tsp_for_dimod_bqm(ncity: int, distances: np.ndarray, dmax: float):
bqm = dimod.BinaryQuadraticModel(ncity * ncity, dimod.BINARY)
# Objective function
for n in range(ncity):
for i in range(ncity):
for j in range(ncity):
bqm.add_quadratic(
n * ncity + i,
((n + 1) % ncity) * ncity + j,
distances[i, j],
)
# Constraint on each row
for n in range(ncity):
left = [(n * ncity + i, 1) for i in range(ncity)]
bqm.add_linear_equality_constraint(left, dmax, -1)
# Constraint on each column
for i in range(ncity):
left = [(n * ncity + i, 1) for n in range(ncity)]
bqm.add_linear_equality_constraint(left, dmax, -1)
return bqm
class BenchTspDimodBQM:
def create_model(self, ncity: int, distances: np.ndarray, dmax: float):
self.model = tsp_for_dimod_bqm(ncity, distances, dmax)
def to_qubo(self):
self.model.to_qubo()
dimod BQM (symbol math)
import dimod
def tsp_for_dimod_bqm_sym(
ncity: int, distances: np.ndarray, dmax: float
) -> dimod.BinaryQuadraticModel:
bqm = dimod.BinaryQuadraticModel(ncity * ncity, dimod.BINARY)
vars = [
[dimod.Binary(f"{n},{i}") for i in range(ncity)] for n in range(ncity)
]
# Objective function
for n in range(ncity):
for i in range(ncity):
for j in range(ncity):
bqm += distances[i, j] * vars[n][i] * vars[(n + 1) % ncity][j]
# Constraint on each row
for n in range(ncity):
bqm += dmax * (sum(vars[n][i] for i in range(ncity)) - 1) ** 2
# Constraint on each column
for i in range(ncity):
bqm += dmax * (sum(vars[n][i] for n in range(ncity)) - 1) ** 2
return bqm # type: ignore
class BenchTspDimodBQMSym:
def create_model(self, ncity: int, distances: np.ndarray, dmax: float):
self.model = tsp_for_dimod_bqm_sym(ncity, distances, dmax)
def to_qubo(self):
self.model.to_qubo()
dimod CQM
import dimod
def tsp_for_dimod_cqm(ncity: int, distances: np.ndarray, dmax: float):
cqm = dimod.ConstrainedQuadraticModel()
vars = [
[dimod.Binary(f"{n},{i}") for i in range(ncity)] for n in range(ncity)
]
# Objective function
obj = 0.0
for n in range(ncity):
for i in range(ncity):
for j in range(ncity):
obj += distances[i, j] * vars[n][i] * vars[(n + 1) % ncity][j]
cqm.set_objective(obj)
# Constraint on each row
for n in range(ncity):
cqm.add_constraint(sum(vars[n]) == 1)
# Constraint on each column
for i in range(ncity):
cqm.add_constraint(sum(vars[n][i] for n in range(ncity)) == 1)
return cqm
class BenchTspDimodCQM:
def create_model(self, ncity: int, distances: np.ndarray, dmax: float):
self.model = tsp_for_dimod_cqm(ncity, distances, dmax)
def to_qubo(self):
pass
Benchmark code
import time
def make_distance(ncity: int) -> tuple[np.ndarray, float]:
rng = np.random.default_rng(12345)
x = rng.random(ncity)
y = rng.random(ncity)
distances = (
(x[:, np.newaxis] - x[np.newaxis, :]) ** 2
+ (y[:, np.newaxis] - y[np.newaxis, :]) ** 2
) ** 0.5
dmax: float = np.max(distances) # type: ignore
return distances, dmax
for ncity in [32, 100, 317]:
distances, dmax = make_distance(ncity)
for bench_class in [
BenchTspAmplify,
BenchTspPyQubo,
BenchTspDimodBQM,
BenchTspDimodBQMSym,
BenchTspDimodCQM,
]:
bench = bench_class()
start = time.time()
bench.create_model(ncity, distances, dmax)
end = time.time()
t1 = end - start
start = time.time()
bench.to_qubo()
end = time.time()
t2 = end - start
print(f"{t1} {t2}")
Benchmark environment
IntelR CoreTM i9-12900K
(E-Cores disabled)
Ubuntu 22.04
amplify 1.0.0
pyqubo 1.4.0
dimod 0.12.14