13from itertools
import product
16from gams
import GamsWorkspace, GamsModifier, SymbolUpdateType
19$if not set tspdata $abort 'tspdata not set'
23 i(ii)
'subset of cities';
25Alias (ii,jj), (i,j,k);
27Parameter c(ii,jj)
'distance matrix';
33$
if not set nrCities $set nrCities 20
34i(ii)$(ord(ii) < %nrCities%) = yes;
37 x(ii,jj)
'decision variables - leg of trip'
38 z
'objective variable';
44 objective
'total cost'
45 rowsum(ii)
'leave each city only once'
46 colsum(jj)
'arrive at each city only once';
48* the assignment problem
is a relaxation of the TSP
49objective.. z =e= sum((i,j), c(i,j)*x(i,j));
51rowsum(i).. sum(j, x(i,j)) =e= 1;
53colsum(j).. sum(i, x(i,j)) =e= 1;
55$
if not set cmax $set cmax 2
59 acut(cut,ii,jj)
'cut constraint matrix'
60 rhscut(cut)
'cut constraint rhs';
62Equation sscut(cut)
'sub_tour elimination cut';
63sscut(cut).. sum((i,j),
acut(cut,i,j)*x(i,j)) =l=
rhscut(cut);
65set
cc(cut)
'previous cuts';
68$
if set cutdata execute_load
'%cutdata%', cc, acut, rhscut;
70acut(cut,i,j)$(
not cc(cut)) = eps;
79if __name__ == "__main__":
80 sys_dir = sys.argv[1] if len(sys.argv) > 1
else None
81 ws = GamsWorkspace(system_directory=sys_dir)
93 cut_data = ws.add_database()
94 cc = cut_data.add_set(
"cc", 1)
95 acut = cut_data.add_parameter(
"acut", 3)
96 rhscut = cut_data.add_parameter(
"rhscut", 1)
102 cmax = cur_cut + cuts_per_round
106 job = ws.add_job_from_string(GAMS_MODEL)
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"] = (
114 + os.path.join(ws.system_directory,
"apifiles",
"Data",
"tsp.gdx")
117 job.run(gams_options=opt, checkpoint=cp)
121 for i
in job.out_db[
"i"]:
125 mi = cp.add_modelinstance()
126 mi_acut = mi.sync_db.add_parameter(
"acut", 3)
127 mi_rhscut = mi.sync_db.add_parameter(
"rhscut", 1)
129 "assign use mip min z", [GamsModifier(mi_acut), GamsModifier(mi_rhscut)]
133 mi.solve(SymbolUpdateType.Accumulate)
134 mi.sync_db[
"acut"].clear()
135 mi.sync_db[
"rhscut"].clear()
139 for i, j
in product(n, n):
140 if mi.sync_db[
"x"].find_record([i, j]).level > 0.5:
144 not_visited = list(n)
148 while graph[i] != not_visited[0]:
151 not_visited = list(filter(
lambda x: x
not in sub_tour, not_visited))
154 for i, j
in product(sub_tour, sub_tour):
155 acut.add_record([f
"c{cur_cut}", i, j]).value = 1
156 mi_acut.add_record([f
"c{cur_cut}", i, j]).value = 1
157 rhscut.add_record(f
"c{cur_cut}").value = len(sub_tour) - 0.5
158 mi_rhscut.add_record([f
"c{cur_cut}"]).value = len(sub_tour) - 0.5
159 cc.add_record(f
"c{cur_cut}")
163 if len(sub_tour) == len(n):
164 print(f
"z={mi.sync_db['z'].first_record().level}")
167 print(f
"{i} -> ", end=
"")