tsp.py
Go to the documentation of this file.
12 
13 from gams import *
14 from itertools import product
15 import os
16 import sys
17 
19  return '''
20 $Title Traveling Salesman Problem Instance with Python
21 
22 $Ontext
23 
24 The sub_tour elimination constraints are generated by a Python
25 script. The MIP is solved over and over, but GAMS have to
26 generate the model only after n cuts have been added.
27 
28 $Offtext
29 
30 $if not set tspdata $abort 'tspdata not set'
31 
32 set ii cities
33  i(ii) subset of cities
34 alias (ii,jj),(i,j,k);
35 
36 parameter c(ii,jj) distance matrix;
37 
38 $gdxin %tspdata%
39 $load ii c
40 
41 $if not set nrCities $set nrCities 20
42 i(ii)$(ord(ii) < %nrCities%) = yes;
43 
44 variables x(ii,jj) decision variables - leg of trip
45  z objective variable;
46 binary variable x; x.fx(ii,ii) = 0;
47 
48 equations objective total cost
49  rowsum(ii) leave each city only once
50  colsum(jj) arrive at each city only once;
51 
52 * the assignment problem is a relaxation of the TSP
53 objective.. z =e= sum((i,j), c(i,j)*x(i,j));
54 rowsum(i).. sum(j, x(i,j)) =e= 1;
55 colsum(j).. sum(i, x(i,j)) =e= 1;
56 
57 $if not set cmax $set cmax 2
58 set cut /c0*c%cmax%/;
59 parameter
60  acut(cut,ii,jj) cut constraint matrix
61  rhscut(cut) cut constraint rhs;
62 
63 equation sscut(cut) sub_tour elimination cut;
64 sscut(cut).. sum((i,j), Acut(cut,i,j)*x(i,j)) =l= RHScut(cut);
65 
66 set cc(cut) previous cuts; cc(cut) = no;
67 $if set cutdata execute_load '%cutdata%', cc, Acut, RHScut;
68 
69 Acut(cut,i,j)$(not cc(cut)) = eps;
70 RHScut(cut)$(not cc(cut)) = card(ii);
71 
72 model assign /all/;
73 
74 option optcr=0; '''
75 
76 
77 if __name__ == "__main__":
78  if len(sys.argv) > 1:
79  ws = GamsWorkspace(system_directory = sys.argv[1])
80  else:
81  ws = GamsWorkspace()
82 
83  # number of cuts that can be added to a model instance
84  cuts_per_round = 10
85  # current cut
86  curcut = 0
87  # cut limit for current model instance (cmax = curcut + cuts_per_round)
88  cmax = 0
89 
90  # database used to collect all generated cuts
91  cut_data = ws.add_database()
92  cc = cut_data.add_set("cc", 1, "")
93  acut = cut_data.add_parameter("acut", 3, "")
94  rhscut = cut_data.add_parameter("rhscut", 1, "")
95  # list of cities (i1, i2, i3, ...)
96  n = []
97 
98  # do-while loop
99  while True:
100  # create a new model instance when the cut limit is reached
101  if curcut >= cmax:
102  cmax = curcut + cuts_per_round
103  cut_data.export();
104 
105  # create the checkpoint
106  tsp_job = ws.add_job_from_string(get_model_text())
107  cp = ws.add_checkpoint()
108  opt = ws.add_options()
109  opt.defines["nrcities"] = "20"
110  opt.defines["cmax"] = str(cmax - 1)
111  opt.defines["cutdata"] = cut_data.name
112  opt.defines["tspdata"] = '"' + os.path.join(os.path.dirname(os.path.realpath(__file__)), "..", "Data", "tsp.gdx") + '"'
113  tsp_job.run(gams_options=opt, checkpoint=cp)
114 
115  # fill the n list only once
116  if not n:
117  for i in tsp_job.out_db["i"]:
118  n += i.keys
119 
120  #instantiate the model instance with modifiers mi_acut and mi_rhscut
121  mi = cp.add_modelinstance()
122  mi_acut = mi.sync_db.add_parameter("acut", 3, "")
123  mi_rhscut = mi.sync_db.add_parameter("rhscut", 1, "")
124  mi.instantiate("assign use mip min z", [GamsModifier(mi_acut), GamsModifier(mi_rhscut)])
125 
126  #solve model instance using update type accumulate and clear acut and rhscut afterwards
127  mi.solve(SymbolUpdateType.Accumulate)
128  mi.sync_db["acut"].clear()
129  mi.sync_db["rhscut"].clear()
130 
131  # collect graph information from the solution
132  graph = {}
133  not_visited = list(n)
134  for i,j in product(n,n):
135  if mi.sync_db["x"].find_record([i, j]).level > 0.5:
136  graph[i] = j
137 
138  # find all subtours and add the required cuts by modifying acut and rhscut
139  while not_visited:
140  i = not_visited[0]
141  sub_tour = [i]
142  while graph[i] != not_visited[0]:
143  i = graph[i]
144  sub_tour.append(i)
145  not_visited = list(filter(lambda x: x not in sub_tour, not_visited))
146 
147  # add the cuts to both databases (cut_data, mi.sync_db)
148  for i,j in product(sub_tour,sub_tour):
149  acut.add_record(["c"+str(curcut), i, j]).value = 1
150  mi_acut.add_record(["c"+str(curcut), i, j]).value = 1
151  rhscut.add_record("c"+str(curcut)).value = len(sub_tour)-0.5
152  mi_rhscut.add_record(["c"+str(curcut)]).value = len(sub_tour)-0.5
153  cc.add_record("c"+str(curcut))
154  curcut += 1
155 
156  # quit the do/while loop if a subtour contains all cities and print the solution
157  if len(sub_tour) == len(n):
158  print("z=" + str(mi.sync_db["z"].first_record().level))
159  print("sub_tour: ")
160  for i in sub_tour:
161  print(i + " -> ", end="")
162  print(sub_tour[0])
163  break;
164 
def get_model_text()
Definition: tsp.py:18