/
tsp.py
241 lines (181 loc) · 7.68 KB
/
tsp.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
# Copyright 2018 D-Wave Systems 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.
import itertools
from collections import defaultdict
from dwave_networkx.utils import binary_quadratic_model_sampler
__all__ = ["traveling_salesperson",
"traveling_salesperson_qubo",
"traveling_salesman",
"traveling_salesman_qubo",
"is_hamiltonian_path",
]
@binary_quadratic_model_sampler(1)
def traveling_salesperson(G, sampler=None, lagrange=None, weight='weight',
start=None, **sampler_args):
"""Returns an approximate minimum traveling salesperson route.
Defines a QUBO with ground states corresponding to the
minimum routes and uses the sampler to sample
from it.
A route is a cycle in the graph that reaches each node exactly once.
A minimum route is a route with the smallest total edge weight.
Parameters
----------
G : NetworkX graph
The graph on which to find a minimum traveling salesperson route.
This should be a complete graph with non-zero weights on every edge.
sampler :
A binary quadratic model sampler. A sampler is a process that
samples from low energy states in models defined by an Ising
equation or a Quadratic Unconstrained Binary Optimization
Problem (QUBO). A sampler is expected to have a 'sample_qubo'
and 'sample_ising' method. A sampler is expected to return an
iterable of samples, in order of increasing energy. If no
sampler is provided, one must be provided using the
`set_default_sampler` function.
lagrange : number, optional (default None)
Lagrange parameter to weight constraints (visit every city once)
versus objective (shortest distance route).
weight : optional (default 'weight')
The name of the edge attribute containing the weight.
start : node, optional
If provided, the route will begin at `start`.
sampler_args :
Additional keyword parameters are passed to the sampler.
Returns
-------
route : list
List of nodes in order to be visited on a route
Examples
--------
>>> import dimod
...
>>> G = nx.Graph()
>>> G.add_weighted_edges_from({(0, 1, .1), (0, 2, .5), (0, 3, .1), (1, 2, .1),
... (1, 3, .5), (2, 3, .1)})
>>> dnx.traveling_salesperson(G, dimod.ExactSolver(), start=0) # doctest: +SKIP
[0, 1, 2, 3]
Notes
-----
Samplers by their nature may not return the optimal solution. This
function does not attempt to confirm the quality of the returned
sample.
"""
# Get a QUBO representation of the problem
Q = traveling_salesperson_qubo(G, lagrange, weight)
# use the sampler to find low energy states
response = sampler.sample_qubo(Q, **sampler_args)
sample = response.first.sample
route = [None]*len(G)
for (city, time), val in sample.items():
if val:
route[time] = city
if start is not None and route[0] != start:
# rotate to put the start in front
idx = route.index(start)
route = route[idx:] + route[:idx]
return route
traveling_salesman = traveling_salesperson
def traveling_salesperson_qubo(G, lagrange=None, weight='weight', missing_edge_weight=None):
"""Return the QUBO with ground states corresponding to a minimum TSP route.
If :math:`|G|` is the number of nodes in the graph, the resulting qubo will have:
* :math:`|G|^2` variables/nodes
* :math:`2 |G|^2 (|G| - 1)` interactions/edges
Parameters
----------
G : NetworkX graph
A complete graph in which each edge has a attribute giving its weight.
lagrange : number, optional (default None)
Lagrange parameter to weight constraints (no edges within set)
versus objective (largest set possible).
weight : optional (default 'weight')
The name of the edge attribute containing the weight.
missing_edge_weight : number, optional (default None)
For bi-directional graphs, the weight given to missing edges.
If None is given (the default), missing edges will be set to
the sum of all weights.
Returns
-------
QUBO : dict
The QUBO with ground states corresponding to a minimum travelling
salesperson route. The QUBO variables are labelled `(c, t)` where `c`
is a node in `G` and `t` is the time index. For instance, if `('a', 0)`
is 1 in the ground state, that means the node 'a' is visted first.
"""
N = G.number_of_nodes()
if lagrange is None:
# If no lagrange parameter provided, set to 'average' tour length.
# Usually a good estimate for a lagrange parameter is between 75-150%
# of the objective function value, so we come up with an estimate for
# tour length and use that.
if G.number_of_edges()>0:
lagrange = G.size(weight=weight)*G.number_of_nodes()/G.number_of_edges()
else:
lagrange = 2
# calculate default missing_edge_weight if required
if missing_edge_weight is None:
# networkx method to calculate sum of all weights
missing_edge_weight = G.size(weight=weight)
# some input checking
if N in (1, 2):
msg = "graph must have at least 3 nodes or be empty"
raise ValueError(msg)
# Creating the QUBO
Q = defaultdict(float)
# Constraint that each row has exactly one 1
for node in G:
for pos_1 in range(N):
Q[((node, pos_1), (node, pos_1))] -= lagrange
for pos_2 in range(pos_1+1, N):
Q[((node, pos_1), (node, pos_2))] += 2.0*lagrange
# Constraint that each col has exactly one 1
for pos in range(N):
for node_1 in G:
Q[((node_1, pos), (node_1, pos))] -= lagrange
for node_2 in set(G)-{node_1}:
# QUBO coefficient is 2*lagrange, but we are placing this value
# above *and* below the diagonal, so we put half in each position.
Q[((node_1, pos), (node_2, pos))] += lagrange
# Objective that minimizes distance
for u, v in itertools.combinations(G.nodes, 2):
for pos in range(N):
nextpos = (pos + 1) % N
# going from u -> v
try:
value = G[u][v][weight]
except KeyError:
value = missing_edge_weight
Q[((u, pos), (v, nextpos))] += value
# going from v -> u
try:
value = G[v][u][weight]
except KeyError:
value = missing_edge_weight
Q[((v, pos), (u, nextpos))] += value
return Q
traveling_salesman_qubo = traveling_salesperson_qubo
def is_hamiltonian_path(G, route):
"""Determines whether the given list forms a valid TSP route.
A travelling salesperson route must visit each city exactly once.
Parameters
----------
G : NetworkX graph
The graph on which to check the route.
route : list
List of nodes in the order that they are visited.
Returns
-------
is_valid : bool
True if route forms a valid travelling salesperson route.
"""
return (set(route) == set(G))