# INTRODUCTION

1. Input the data and visualize the problem
2. Model TSP as ILP formulation w/o Subtour constraints
3. Implement Subtour Elimination Method 1: MTZ’s Approach
4. Implement Subtour Elimination Method 2: DFJ’s Approach
5. Compare MTZ’s formulation and DFJ’s formulation
6. Conclusion

# 1 Input the data and problem visualization

`#import libraries%matplotlib inlineimport pulpimport pandas as pdfrom scipy.spatial import distance_matrixfrom matplotlib import pyplot as pltimport timeimport copy`
`#This function takes locations as input and plot a scatter plotdef plot_fig(loc,heading="plot"):plt.figure(figsize=(10,10))    for i,row in loc.iterrows():        if i==0:            plt.scatter(row["x"],row["y"],c='r')            plt.text(row["x"]+0.2, row["y"]+0.2, 'DELHI (depot) ')        else:            plt.scatter(row["x"], row["y"], c='black')            plt.text(row["x"] + 0.2, row["y"] + 0.2,full_data.loc[i]['CITy'] )        plt.ylim(6,36)        plt.xlim(66,96)        plt.title(heading)# this function find all the subtour in the LP solution.def get_plan(r0):    r=copy.copy(r0)    route = []    while len(r) != 0:        plan = [r]        del (r)        l = 0        while len(plan) > l:            l = len(plan)            for i, j in enumerate(r):                if plan[-1] == j:                    plan.append(j)                    del (r[i])        route.append(plan)    return(route)`
`# set no of citiesno_of_locs=6data=pd.read_csv("tsp_city_data.csv")full_data=data.iloc[0:no_of_locs,:]d=full_data[['x','y']]dis_mat=pd.DataFrame(distance_matrix(d.values,d.values),\                       index=d.index,columns=d.index)print("----------data--------------")print(full_data)print("-----------------------------")plot_fig(d,heading="Problem Visualization")   plt.show()`

# 2 Model TSP in ILP without Subtour elimination constraints

`model=pulp.LpProblem('tsp',pulp.LpMinimize)#define variablex=pulp.LpVariable.dicts("x",((i,j) for i in range(no_of_locs) \                             for j in range(no_of_locs)),\                            cat='Binary')#set objectivemodel+=pulp.lpSum(dis_mat[i][j]* x[i,j] for i in range(no_of_locs) \                      for j in range(no_of_locs))# st constraintsfor i in range(no_of_locs):    model+=x[i,i]==0    model+=pulp.lpSum(x[i,j] for j in range(no_of_locs))==1    model += pulp.lpSum(x[j, i] for j in range(no_of_locs)) == 1status=model.solve()#status=model.solver()print("-----------------")print(status,pulp.LpStatus[status],pulp.value(model.objective))route=[(i,j) for i in range(no_of_locs) \           for j in range(no_of_locs) if pulp.value(x[i,j])==1]print(get_plan(route))`
`plot_fig(d,heading="solution Visualization")arrowprops = dict(arrowstyle='->', connectionstyle='arc3', edgecolor='blue')for i, j in route:    plt.annotate('', xy=[d.iloc[j]['x'], d.iloc[j]['y']],\                      xytext=[d.iloc[i]['x'], d.iloc[i]['y']],\                     arrowprops=arrowprops)`
1. Tour 1 : Delhi > Nagpur > Rajkot
2. Tour 2 : Kolkata > Dispur > Agartala

# How does constraint (5) remove subtours ?

`start_t=time.time()model=pulp.LpProblem('tsp',pulp.LpMinimize)#define variablex=pulp.LpVariable.dicts("x",((i,j) for i in range(no_of_locs) \                                for j in range(no_of_locs)),\                           cat='Binary')t = pulp.LpVariable.dicts("t", (i for i in range(no_of_locs)), \                             lowBound=1,upBound= no_of_locs, cat='Continuous')#set objectivemodel+=pulp.lpSum(dis_mat[i][j]* x[i,j] for i in range(no_of_locs) \                      for j in range(no_of_locs))# st constraintsfor i in range(no_of_locs):    model+=x[i,i]==0    model+=pulp.lpSum(x[i,j] for j in range(no_of_locs))==1    model += pulp.lpSum(x[j, i] for j in range(no_of_locs)) == 1#eliminate subtourfor i in range(no_of_locs):    for j in range(no_of_locs):        if i!=j and (i!=0 and j!=0):            model+=t[j]>=t[i]+1 - (2*no_of_locs)*(1-x[i,j])status=model.solve()#status=model.solver()print("-----------------")print(status,pulp.LpStatus[status],pulp.value(model.objective))route=[(i,j) for i in range(no_of_locs) \           for j in range(no_of_locs) if pulp.value(x[i,j])==1]print(route)plot_fig(d,heading="solution Visualization")arrowprops = dict(arrowstyle='->', connectionstyle='arc3', edgecolor='blue')for i, j in route:    plt.annotate('', xy=[d.iloc[j]['x'], d.iloc[j]['y']],\                      xytext=[d.iloc[i]['x'], d.iloc[i]['y']],\                     arrowprops=arrowprops)print("time taken by MTZ formulation = ", time.time()-start_t)`

# Higher-level Algorithm for DFJ

`start_t_1=time.time()model=pulp.LpProblem('tsp',pulp.LpMinimize)#define variablex=pulp.LpVariable.dicts("x",((i,j) for i in range(no_of_locs) \                                for j in range(no_of_locs)),\                           cat='Binary')#set objectivemodel+=pulp.lpSum(dis_mat[i][j]* x[i,j] for i in range(no_of_locs) \                      for j in range(no_of_locs))# st constraintsfor i in range(no_of_locs):    model+=x[i,i]==0    model+=pulp.lpSum(x[i,j] for j in range(no_of_locs))==1    model += pulp.lpSum(x[j, i] for j in range(no_of_locs)) == 1    status=model.solve()route=[(i,j) for i in range(no_of_locs) \           for j in range(no_of_locs) if pulp.value(x[i,j])==1]route_plan=get_plan(route)subtour=[]while len(route_plan)!=1:    for i in range(len(route_plan)):            #print(route_plan[i])model+=pulp.lpSum(x[route_plan[i][j],route_plan[i][j]]\                              for j in range(len(route_plan[i])))<=\                              len(route_plan[i])-1    status=model.solve()    route = [(i, j) for i in range(no_of_locs) \                 for j in range(no_of_locs) if pulp.value(x[i, j]) == 1]route_plan = get_plan(route)        subtour.append(len(route_plan))print("-----------------")print(status,pulp.LpStatus[status],pulp.value(model.objective))print(route_plan)print("no. of times LP model is solved = ",len(subtour))print("subtour log (no. of subtours in each solution))",subtour)print("Time taken by DFJ formulation = ", time.time()-start_t_1)plot_fig(d,heading="solution Visualization")arrowprops = dict(arrowstyle='->', connectionstyle='arc3', edgecolor='blue')for i, j in route_plan:    plt.annotate('', xy=[d.iloc[j]['x'], d.iloc[j]['y']],\                     xytext=[d.iloc[i]['x'], d.iloc[i]['y']],                     arrowprops=arrowprops)plt.show()#print("total time = ",time.time()-start)`

# Conclusion

--

--

## More from The Startup

Get smarter at building your thing. Follow to join The Startup’s +8 million monthly readers & +768K followers.