Repository: dwavesystems/dimod
Pull request #1380 (2024-07-19 22:51:53)

create qap generator


Patch:
diff --git a/dimod/generators/__init__.py b/dimod/generators/__init__.py
index 5db1b293..d6e9ac06 100644
--- a/dimod/generators/__init__.py
+++ b/dimod/generators/__init__.py
@@ -24,6 +24,7 @@ from dimod.generators.integer import *
 from dimod.generators.knapsack import *
 from dimod.generators.magic_square import *
 from dimod.generators.multi_knapsack import *
+from dimod.generators.quadratic_assignment import *
 from dimod.generators.random import *
 from dimod.generators.satisfiability import *
 from dimod.generators.wireless import *
\ No newline at end of file
diff --git a/dimod/generators/quadratic_assignment.py b/dimod/generators/quadratic_assignment.py
new file mode 100644
index 00000000..37143f39
--- /dev/null
+++ b/dimod/generators/quadratic_assignment.py
@@ -0,0 +1,79 @@
+# Copyright 2024 D-Wave Inc.
+#
+#    Licensed under the Apache License, Version 2.0 (the "License");
+#    you may not use this file except in compliance with the License.
+#    You may obtain a copy of the License at
+#
+#        http://www.apache.org/licenses/LICENSE-2.0
+#
+#    Unless required by applicable law or agreed to in writing, software
+#    distributed under the License is distributed on an "AS IS" BASIS,
+#    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+#    See the License for the specific language governing permissions and
+#    limitations under the License.
+
+from itertools import product
+
+import numpy as np
+
+from dimod.binary_quadratic_model import BinaryQuadraticModel
+from dimod.constrained import ConstrainedQuadraticModel
+from dimod.typing import ArrayLike
+
+__all__ = ['quadratic_assignment']
+
+def quadratic_assignment(distance_matrix: ArrayLike,
+                         flow_matrix: ArrayLike,
+                         ) -> ConstrainedQuadraticModel:
+    """Generates a constrained quadratic model encoding a quadratic-assignment problem.
+
+    Given distance and flow matrices, generates a :class:`~dimod.ConstrainedQuadraticModel` 
+    for the corresponding quadratic-assignment problem. 
+
+    Args:
+        distances: Distances between locations :math:`i` and :math:`j` as a NumPy array.
+
+        flows: Flows between facilities :math:`i` and :math:`j` as a NumPy array.
+
+    Returns:
+        The constrained quadratic model encoding the quadratic-assignment problem. Variables are
+        denoted as ``x_{i}_{j}`` where ``x_{i}_{j} == 1`` means that facility ``i`` is
+        placed in location ``j``.
+
+    """
+
+    distance_matrix = np.atleast_2d(np.asarray(distance_matrix))
+    flow_matrix = np.atleast_2d(np.asarray(flow_matrix))
+
+    if distance_matrix.shape != flow_matrix.shape:
+        raise ValueError("'distance_matrix' and 'flow_matrix' must have the same shape")
+
+    if distance_matrix.shape[0] != distance_matrix.shape[1]:
+        raise ValueError("'distance_matrix' must be square")
+
+    if distance_matrix.ndim != 2:
+        raise ValueError("'distance_matrix' must be 2-dimensional")
+
+    num_locations = distance_matrix.shape[0]
+
+    model = ConstrainedQuadraticModel()
+    obj = BinaryQuadraticModel(vartype='BINARY')
+    x = {(i, j): obj.add_variable(f'x_{i}_{j}') for i in range(num_locations)
+         for j in range(num_locations)}
+
+    for i, j, k, l in product(range(num_locations), repeat=4):
+        if (i,j) != (k,l):
+            obj.set_quadratic(x[(i,j)], x[(k,l)], flow_matrix[i][k]*distance_matrix[j][l]
+                              + flow_matrix[k][i]*distance_matrix[j][l])
+
+    model.set_objective(obj)
+
+    for i in range(num_locations):
+        constraint_vars = [x[(i,j)] for j in range(num_locations)]
+        model.add_discrete(constraint_vars, label=f'discrete_constraint_{i}')
+
+    for j in range(num_locations):
+        constraint = [(x[(i,j)], 1) for i in range(num_locations)] + [(-1,)]
+        model.add_constraint(constraint, sense='==', label=f'facility_constraint_{j}')
+
+    return model
diff --git a/docs/reference/generators.rst b/docs/reference/generators.rst
index c69601d2..0400690a 100644
--- a/docs/reference/generators.rst
+++ b/docs/reference/generators.rst
@@ -48,6 +48,7 @@ Optimization
    maximum_independent_set
    maximum_weight_independent_set
    mimo
+   quadratic_assignment
    random_bin_packing
    random_knapsack
    random_multi_knapsack
diff --git a/releasenotes/notes/qap-generator-0687ed50c85670cd.yaml b/releasenotes/notes/qap-generator-0687ed50c85670cd.yaml
new file mode 100644
index 00000000..fba99e7e
--- /dev/null
+++ b/releasenotes/notes/qap-generator-0687ed50c85670cd.yaml
@@ -0,0 +1,4 @@
+---
+features:
+  - |
+    Add ``quadratic_assignment`` generator function to ``dimod.generators``.
\ No newline at end of file


Test patch:
diff --git a/tests/test_generators.py b/tests/test_generators.py
index c7dcc9db..ec0af0a4 100644
--- a/tests/test_generators.py
+++ b/tests/test_generators.py
@@ -1734,3 +1734,74 @@ class TestPowerR(unittest.TestCase):
         self.assertTrue(all(val != 0 for val in bqm.quadratic.values()))
         self.assertTrue(all(val <= 5 for val in bqm.quadratic.values()))
         self.assertTrue(all(val >= -5 for val in bqm.quadratic.values()))
+
+
+class TestQuadraticAssignment(unittest.TestCase):
+
+    def test_unequal_shape(self):
+        distance_matrix = np.array([[1,2], [3,4]])
+        flow_matrix = np.array([5])
+        with self.assertRaises(ValueError):
+            dimod.generators.quadratic_assignment(distance_matrix, flow_matrix)
+
+    def test_not_square(self):
+        distance_matrix = np.array([[1,2],[3,4],[5,6]])
+        flow_matrix = np.array([[1,2],[3,4],[5,6]])
+        with self.assertRaises(ValueError):
+            dimod.generators.quadratic_assignment(distance_matrix, flow_matrix)
+
+    def test_not_2d(self):
+        distance_matrix = np.array([
+                                    [[1, 2, 3, 4],
+                                    [5, 6, 7, 8]],                                    
+                                    [[13, 14, 15, 16],
+                                    [17, 18, 19, 20]]
+                                ])
+        flow_matrix = np.array([
+                                    [[1, 2, 3, 4],
+                                    [5, 6, 7, 8]],                                    
+                                    [[13, 14, 15, 16],
+                                    [17, 18, 19, 20]]
+                                ])
+        with self.assertRaises(ValueError):
+            dimod.generators.quadratic_assignment(distance_matrix, flow_matrix)
+
+    def test_model(self):
+        distance_matrix = np.array([[0,2,3],[4,0,5],[6,1,0]])
+        flow_matrix = np.array([[0,1,7],[2,0,3],[8,5,0]])
+        num_locations = distance_matrix.shape[0]
+        cqm = dimod.generators.quadratic_assignment(distance_matrix, flow_matrix)
+        self.assertEqual(len(cqm.variables), num_locations**2)
+        self.assertEqual(len(cqm.constraints), 2*num_locations)
+
+    def test_infeasible(self):
+        distance_matrix = np.array([[0,2,3],[4,0,5],[6,1,0]])
+        flow_matrix = np.array([[0,1,7],[2,0,3],[8,5,0]])
+        num_locations = distance_matrix.shape[0]
+        cqm = dimod.generators.quadratic_assignment(distance_matrix, flow_matrix)
+
+        for i in range(num_locations):
+            x = {f'x_{i}_{j}': 1 for j in range(num_locations)}
+            lhs = cqm.constraints[f'discrete_constraint_{i}'].lhs.energy(x)
+            self.assertNotEqual(lhs, 1)
+
+        for j in range(num_locations):
+            x = {f'x_{i}_{j}': 1 for i in range(num_locations)}
+            lhs = cqm.constraints[f'facility_constraint_{j}'].lhs.energy(x)
+            self.assertNotEqual(lhs, 0)
+
+    def test_feasible(self):
+        distance_matrix = np.array([[0,2,3],[4,0,5],[6,1,0]])
+        flow_matrix = np.array([[0,1,7],[2,0,3],[8,5,0]])
+        num_locations = distance_matrix.shape[0]
+        cqm = dimod.generators.quadratic_assignment(distance_matrix, flow_matrix)
+
+        for i in range(num_locations):
+            x = {f'x_{i}_{j}': 1 if i==j else 0 for j in range(num_locations)}
+            lhs = cqm.constraints[f'discrete_constraint_{i}'].lhs.energy(x)
+            self.assertEqual(lhs, 1)
+
+        for j in range(num_locations):
+            x = {f'x_{i}_{j}': 1 if i==j else 0 for i in range(num_locations)}
+            lhs = cqm.constraints[f'facility_constraint_{j}'].lhs.energy(x)
+            self.assertEqual(lhs, 0)