Setup environment

importpandasaspdimportnumpyasnp%matplotlibinline%configInlineBackend.figure_format='retina'importmatplotlib.pyplotaspltfrompylabimportrcParams

/Users/rogerio/.pyenv/versions/3.6.5/lib/python3.6/importlib/_bootstrap.py:219: RuntimeWarning: numpy.dtype size changed, may indicate binary incompatibility. Expected 96, got 88  return f(*args, **kwds)/Users/rogerio/.pyenv/versions/3.6.5/lib/python3.6/importlib/_bootstrap.py:219: RuntimeWarning: numpy.dtype size changed, may indicate binary incompatibility. Expected 96, got 88  return f(*args, **kwds)

Customize graphics style

print(plt.style.available)

['seaborn-dark', 'seaborn-darkgrid', 'seaborn-ticks', 'fivethirtyeight', 'seaborn-whitegrid', 'classic', '_classic_test', 'fast', 'seaborn-talk', 'seaborn-dark-palette', 'seaborn-bright', 'seaborn-pastel', 'grayscale', 'seaborn-notebook', 'ggplot', 'seaborn-colorblind', 'seaborn-muted', 'seaborn', 'Solarize_Light2', 'seaborn-paper', 'bmh', 'tableau-colorblind10', 'seaborn-white', 'dark_background', 'seaborn-poster', 'seaborn-deep']

plt.style.use('classic')rcParams['figure.figsize']=16,10

Initialize data

np.random.seed(18)size,cmin,cmax=150,-100,100data=pd.DataFrame(dict(x=[0],y=[0])).append(pd.DataFrame((np.random.random_sample(2*size)*(cmax-cmin)+cmin).reshape(-1,2),columns=['x','y']),ignore_index=True)data

.dataframe tbody tr th {    vertical-align: top;}.dataframe thead th {    text-align: right;}

151 rows × 2 columns

Plot problem

plt.scatter(data.x,data.y,c='r',marker='o',s=40)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Compute distance matrix

distance_func=lambdaa,b:np.sqrt((a.x-b.x)**2+ (a.y-b.y)**2)defcompute_distances(dfunc=distance_func):foriinrange(size+1):current=data.iloc[i]data[i]=dfunc(current,data)#data[str(i)][i]=np.nan%timecompute_distances()data

CPU times: user 346 ms, sys: 5.93 ms, total: 352 msWall time: 357 ms

.dataframe tbody tr th {    vertical-align: top;}.dataframe thead th {    text-align: right;}

151 rows × 153 columns

Helper functions

estimate_cost=lambdaroute:sum(data.iloc[i][j]fori,jinzip(route,route[1:]))get_coords=lambdaroute:list(zip(*[(data.iloc[i].x,data.iloc[i].y)foriinroute]))

Scenario: center depot, one van, nearest neighbor heuristic

Apply a greedy NN algorithm

defnearest_neighbor(unserved,stop_condition=lambdaroute,target:False):current=0#depotresult_path=[]whileTrue:result_path.append(current)unserved.remove(current)ifnotunserved:breakcurrent=data.iloc[unserved,current+2].idxmin()ifstop_condition(result_path,int(current)):iflen(result_path)>1:breakresult_path.append(0)returnresult_path%timeroute=nearest_neighbor(list(range(size+1)))print('cost={}\nroute={}'.format(estimate_cost(route),route))

CPU times: user 116 ms, sys: 8.29 ms, total: 125 msWall time: 122 mscost=2623.828598616268route=[0, 53, 60, 20, 66, 148, 55, 51, 1, 120, 104, 123, 74, 93, 130, 42, 29, 19, 91, 107, 81, 11, 134, 117, 57, 30, 34, 67, 75, 46, 39, 131, 16, 114, 72, 96, 38, 22, 122, 92, 36, 126, 113, 7, 58, 111, 28, 108, 63, 2, 137, 50, 128, 69, 141, 18, 146, 37, 54, 95, 24, 121, 35, 83, 3, 49, 12, 73, 105, 129, 90, 144, 9, 86, 132, 31, 70, 127, 149, 6, 71, 41, 140, 76, 102, 32, 100, 97, 139, 88, 78, 26, 17, 135, 27, 13, 110, 89, 59, 119, 40, 15, 56, 145, 10, 106, 82, 94, 115, 87, 98, 44, 64, 52, 14, 33, 136, 124, 65, 142, 150, 147, 143, 109, 112, 23, 47, 101, 45, 138, 99, 133, 8, 48, 62, 80, 103, 43, 68, 61, 125, 21, 25, 118, 5, 79, 116, 77, 85, 84, 4, 0]

Plot result path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve it with a 2-opt heuristic

deftwoOptSwap(route,i,j):t=route[i:j+1];t.reverse();route[i:j+1]=tdeftwoOpt(route):cost=estimate_cost(route)definner():nonlocalcostforiinrange(1,len(route)-2):nodei=data.iloc[route[i]]forjinrange(i+1,len(route)-1):nodej=data.iloc[route[j]]save=nodei.iloc[route[i-1]+2]+nodej.iloc[route[j+1]+2]- (nodej.iloc[route[i-1]+2]+nodei.iloc[route[j+1]+2])ifsave>0:twoOptSwap(route,i,j)cost-=saveprint('exchanging {}-{},{}-{} with {}-{},{}-{} => save={} => cost: {}'.format(route[i-1],route[i],route[j],route[j+1],route[i-1],route[j],route[i],route[j+1],save,cost))returnFalsereturnTruewhileTrue:ifinner():break%timetwoOpt(route)

exchanging 53-64,60-52 with 53-60,64-52 => save=2.8803537373675177 => cost: 2620.948244878901...exchanging 1-120,81-107 with 1-81,120-107 => save=26.496627104257605 => cost: 2100.7350743873462CPU times: user 3min 11s, sys: 1.65 s, total: 3min 13sWall time: 3min 21s

Plot improved path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Incredible! A clean path, without any crossing parts!

Scenario: center depot, one van, 2-opt heuristic

Just for fun, let's pass a random route to 2-opt heuristic

route= [xforxinrange(size+1)]+[0]print('cost={}\nroute={}'.format(estimate_cost(route),route))

cost=16332.039263088578route=[0, 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, 0]

Plot result path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's see what 2-opt heuristic can do

%timetwoOpt(route)

exchanging 0-4,1-5 with 0-1,4-5 => save=11.966404491447804 => cost: 16320.07285859713...exchanging 47-23,112-8 with 47-112,23-8 => save=11.556726776627453 => cost: 2181.5435075388827CPU times: user 19min 11s, sys: 8.18 s, total: 19min 19sWall time: 19min 44s

Plot improved path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

WOW, very cool, a clean path too, almost as good as the previous one.

But extremely slower... It took ~19 minutes, while NN + 2-opt took only ~3 minutes.

Scenario: center depot, several vans, no clustering

Let's apply a sequential NN algorithm, but including a max cost constraint

defsequential_nn(max_cost):unserved=list(range(size+1))result=[]whileTrue:result_path=nearest_neighbor(unserved,lambdaroute,target:estimate_cost(route+[target,0])>=max_cost)result.append(result_path)print('van {}: {}'.format(len(result),estimate_cost(result_path)))ifnotunserved:breakunserved.append(0)returnresult%timeresult_paths=sequential_nn(400)

van 1: 389.0742611200774van 2: 398.3533734081084van 3: 397.4195886051539van 4: 390.9917686566429van 5: 371.6762623788105van 6: 320.29683749922486van 7: 364.8394134778564van 8: 320.99556479051967van 9: 210.1189842618746van 10: 276.59421133431306CPU times: user 1.48 s, sys: 24.4 ms, total: 1.51 sWall time: 1.54 s

Plot result paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve them with a 2-opt heuristic

defimprove_seq_nn():fori,pathinenumerate(result_paths):print('Starting path {}/{}'.format(i+1,len(result_paths)))twoOpt(path)%timeimprove_seq_nn()

Starting path 1/10exchanging 55-19,51-91 with 55-51,19-91 => save=2.6873147586143844 => cost: 386.386946361463...exchanging 67-39,75-0 with 67-75,39-0 => save=4.275672113835263 => cost: 368.50842741951476Starting path 2/10exchanging 17-97,78-0 with 17-78,97-0 => save=0.09782392307953103 => cost: 398.25554948502884...exchanging 88-78,26-0 with 88-26,78-0 => save=1.6805548383893267 => cost: 394.249851279636Starting path 3/10Starting path 4/10exchanging 124-127,65-0 with 124-65,127-0 => save=0.5258127779938491 => cost: 390.465955878649exchanging 14-65,33-0 with 14-33,65-0 => save=3.343148155406432 => cost: 387.1228077232426Starting path 5/10exchanging 0-18,141-69 with 0-141,18-69 => save=6.7520080457503155 => cost: 364.9242543330602exchanging 50-121,146-0 with 50-146,121-0 => save=41.908681513869396 => cost: 323.0155728191908Starting path 6/10exchanging 47-45,101-138 with 47-101,45-138 => save=1.4275147235708587 => cost: 318.869322775654Starting path 7/10exchanging 3-35,83-77 with 3-83,35-77 => save=3.633104522812438 => cost: 361.20630895504394...exchanging 90-129,144-0 with 90-144,129-0 => save=5.489882024385494 => cost: 347.690672350374Starting path 8/10Starting path 9/10Starting path 10/10CPU times: user 3.25 s, sys: 37.7 ms, total: 3.29 sWall time: 3.31 s

Plot improved paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Great! Several clean paths, but with many superpositions, let's improve them once more.

Scenario: center depot, several vans, KMeans clustering

Compute delivery clusters

fromsklearn.clusterimportKMeanswork_per_van=20#how many deliveries each van will start, ie granularity; a better approach would be capacity.model=KMeans(n_clusters=1+size//work_per_van,init='k-means++',random_state=0)%timemodel.fit(data[['x','y']])

CPU times: user 75.7 ms, sys: 4.05 ms, total: 79.7 msWall time: 92.8 msKMeans(algorithm='auto', copy_x=True, init='k-means++', max_iter=300,    n_clusters=8, n_init=10, n_jobs=1, precompute_distances='auto',    random_state=0, tol=0.0001, verbose=0)

Plot clusters

plt.scatter(data.x,data.y,marker='o',c=model.labels_,s=200)plt.scatter(model.cluster_centers_[:,0],model.cluster_centers_[:,1],marker='+',c='r',s=500)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's apply a greedy NN algorithm sequentially, on each of the clusters

defgreedy_nn_clusters():result=[]foriinrange(model.n_clusters):unserved=data[model.labels_==i].index.tolist()if0notinunserved:unserved.append(0)result_path=nearest_neighbor(unserved)result.append(result_path)print('van {}: {}'.format(i+1,estimate_cost(result_path)))returnresult%timeresult_paths=greedy_nn_clusters()

van 1: 409.3080162039939van 2: 266.4745473994821van 3: 413.3606064604762van 4: 446.80922078430706van 5: 457.66477785115126van 6: 455.80304195030294van 7: 459.76894476196577van 8: 376.5008949566861CPU times: user 164 ms, sys: 10.1 ms, total: 174 msWall time: 174 ms

Plot result paths

plt.scatter(data.x,data.y,marker='o',c=model.labels_,s=200)plt.scatter(model.cluster_centers_[:,0],model.cluster_centers_[:,1],marker='+',c='r',s=500)forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],'k-',linewidth=2)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve them with a 2-opt heuristic

defimprove_result_paths():fori,pathinenumerate(result_paths):print('Starting path {}/{}'.format(i+1,len(result_paths)))twoOpt(path)%timeimprove_result_paths()

Starting path 1/8exchanging 48-112,99-138 with 48-99,112-138 => save=11.5162692633929 => cost: 397.791746940601exchanging 45-118,25-5 with 45-25,118-5 => save=13.984640981159586 => cost: 383.8071059594414Starting path 2/8exchanging 148-52,136-65 with 148-136,52-65 => save=1.6057804518926915 => cost: 264.8687669475894exchanging 148-124,52-142 with 148-52,124-142 => save=6.904413311483367 => cost: 257.964353636106Starting path 3/8exchanging 63-24,108-111 with 63-108,24-111 => save=11.657922849361384 => cost: 401.7026836111148...exchanging 113-36,126-0 with 113-126,36-0 => save=7.772149456496223 => cost: 378.1902685773841Starting path 4/8exchanging 147-109,150-0 with 147-150,109-0 => save=5.601665728075147 => cost: 441.2075550562319...exchanging 62-109,143-150 with 62-143,109-150 => save=0.5147852919035643 => cost: 394.93412566866925Starting path 5/8exchanging 9-149,6-144 with 9-6,149-144 => save=0.0011015802705607314 => cost: 457.66367627088067...exchanging 105-129,144-90 with 105-144,129-90 => save=3.918490824545831 => cost: 407.2051777421843Starting path 6/8exchanging 30-69,34-116 with 30-34,69-116 => save=1.1688372549307786 => cost: 454.63420469537215...exchanging 22-116,79-38 with 22-79,116-38 => save=5.039989423215701 => cost: 380.39750649818404Starting path 7/8exchanging 0-55,1-107 with 0-1,55-107 => save=14.962441964752642 => cost: 444.8065027972131...exchanging 0-107,55-42 with 0-55,107-42 => save=6.044195398723886 => cost: 403.197665386089Starting path 8/8exchanging 0-71,26-59 with 0-26,71-59 => save=17.738760993960298 => cost: 358.76213396272584...exchanging 17-78,26-88 with 17-26,78-88 => save=0.17155006280330554 => cost: 319.53217682612376CPU times: user 2.81 s, sys: 64.8 ms, total: 2.87 sWall time: 3.03 s

Plot improved paths

plt.scatter(data.x,data.y,marker='o',c=model.labels_,s=200)plt.scatter(model.cluster_centers_[:,0],model.cluster_centers_[:,1],marker='+',c='r',s=500)forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],'k-',linewidth=2)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Cool clean paths, with clearly defined cluster subsections!

Let's change the depot to the farthest point

Find farthest point and set it as new depot

maxid=data[0].idxmax()print('farthest point={}; distance={}; coords={}'.format(maxid,data.iloc[maxid,2],tuple(data.iloc[maxid][['x','y']].tolist())))data.drop(data.columns.difference(['x','y']),axis=1,inplace=True)data.iloc[0]=data.iloc[maxid]data.drop(maxid,inplace=True)data.reset_index(drop=True,inplace=True)size-=1data

farthest point=84; distance=138.29710566715653; coords=(99.07656719022526, 96.48794364952255)

.dataframe tbody tr th {    vertical-align: top;}.dataframe thead th {    text-align: right;}

150 rows × 2 columns

Plot relocated depot

plt.scatter(data.x,data.y,c='r',marker='o',s=40)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Recompute distance matrix

%timecompute_distances()data

CPU times: user 365 ms, sys: 5.95 ms, total: 371 msWall time: 374 ms

.dataframe tbody tr th {    vertical-align: top;}.dataframe thead th {    text-align: right;}

150 rows × 152 columns

Scenario: farthest depot, one van, nearest neighbor heuristic

Apply a greedy NN algorithm

%timeroute=nearest_neighbor(list(range(size+1)))print('cost={}\nroute={}'.format(estimate_cost(route),route))

CPU times: user 111 ms, sys: 13.4 ms, total: 124 msWall time: 118 mscost=2598.9001149206997route=[0, 84, 77, 83, 3, 49, 12, 73, 104, 128, 89, 143, 9, 85, 131, 31, 70, 126, 148, 6, 71, 41, 139, 76, 101, 32, 99, 96, 138, 87, 78, 26, 17, 134, 27, 13, 109, 88, 59, 118, 40, 15, 56, 144, 10, 105, 82, 93, 114, 86, 97, 44, 64, 53, 60, 20, 66, 147, 55, 51, 1, 119, 103, 122, 74, 92, 129, 42, 29, 19, 90, 106, 81, 11, 133, 116, 57, 30, 34, 67, 75, 46, 39, 130, 16, 113, 72, 95, 38, 22, 121, 91, 36, 125, 112, 7, 58, 110, 28, 107, 63, 2, 136, 50, 127, 69, 140, 18, 145, 37, 54, 94, 24, 120, 35, 123, 135, 33, 14, 52, 65, 141, 149, 146, 142, 108, 111, 23, 47, 100, 45, 137, 98, 132, 8, 48, 62, 80, 102, 43, 68, 61, 124, 21, 25, 117, 5, 79, 115, 4, 0]

Plot result path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve it with a 2-opt heuristic

%timetwoOpt(route)

exchanging 84-4,77-0 with 84-77,4-0 => save=15.967192548480725 => cost: 2582.932922372219...exchanging 3-49,12-83 with 3-12,49-83 => save=1.0654890882007315 => cost: 2111.9918767877866CPU times: user 3min 36s, sys: 1.67 s, total: 3min 38sWall time: 3min 44s

Plot improved path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

A great optimized route too.

Scenario: farthest depot, one van, 2-opt heuristic

Just for fun, let's pass a random route to 2-opt heuristic

route= [xforxinrange(size+1)]+[0]print('cost={}\nroute={}'.format(estimate_cost(route),route))

cost=16502.514122712284route=[0, 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, 0]

Plot result path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's see what 2-opt heuristic can do

%timetwoOpt(route)

exchanging 0-2,1-3 with 0-1,2-3 => save=6.263294946633721 => cost: 16496.25082776565...exchanging 31-70,126-131 with 31-126,70-131 => save=1.7915970025161734 => cost: 2147.6785119636334CPU times: user 30min 1s, sys: 27.7 s, total: 30min 29sWall time: 33min 24s

Plot improved path

coords=get_coords(route)plt.plot(coords[0],coords[1],'b-',linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Yes, it's beautiful what this heuristic can do.

Only again extremely slower. It's best suited to improve a path, not to build one.

Scenario: farthest depot, several vans, no clustering

Let's apply a sequential NN algorithm, but including a max cost constraint

%timeresult_paths=sequential_nn(400)

van 1: 385.64822798111754van 2: 396.65132151762236van 3: 356.93203827425646van 4: 396.0161561312648van 5: 387.061496278307van 6: 395.6147269472501van 7: 369.199156243418van 8: 394.6734767469922van 9: 397.84112119751114van 10: 398.459895479128van 11: 398.88707002367426van 12: 393.1275886242274van 13: 399.0252974546198van 14: 369.5354695202354van 15: 399.04095929693733van 16: 398.1246462968952van 17: 390.2323088432243van 18: 391.6922432861831van 19: 396.62966713208937van 20: 397.37088575736055van 21: 398.69882895151926van 22: 398.71167524672734van 23: 401.9265755043099van 24: 405.5614801848182van 25: 405.6059533951813van 26: 406.6405588064364van 27: 409.6579507334966van 28: 416.22815658152126van 29: 426.0013563494297van 30: 427.5654816187484van 31: 428.5781354316992van 32: 438.7764370279704van 33: 439.2553078706148van 34: 445.7400163235422van 35: 447.1961602342642van 36: 461.09599384348525van 37: 464.4374453992292van 38: 471.7901899743616van 39: 482.2725996546565van 40: 482.59749534967307van 41: 483.3659132930797van 42: 494.26350457825447van 43: 525.5364809069648van 44: 551.3280289580272CPU times: user 1.22 s, sys: 48.8 ms, total: 1.27 sWall time: 1.43 s

Plot result paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Max cost 400 doesn't quite work here, as the distances are much larger

%timeresult_paths=sequential_nn(500)

van 1: 493.1613509811791van 2: 498.37402590912467van 3: 499.57415652250484van 4: 490.0564767221547van 5: 497.67366460125766van 6: 491.0474410501304van 7: 476.1857519522164van 8: 487.87440569372217van 9: 462.17325815341474van 10: 486.2976818275437van 11: 498.3091236408933van 12: 466.4397450832517van 13: 464.4374453992292van 14: 494.26350457825447van 15: 525.5364809069648van 16: 551.3280289580272CPU times: user 573 ms, sys: 20.1 ms, total: 593 msWall time: 588 ms

Plot result paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve them with a 2-opt heuristic

%timeimprove_result_paths()

Starting path 1/16exchanging 84-134,77-0 with 84-77,134-0 => save=17.58671573136749 => cost: 475.57463524981165...exchanging 131-9,85-73 with 131-85,9-73 => save=0.22871091978264957 => cost: 415.7579280346732Starting path 2/16exchanging 114-10,93-0 with 114-93,10-0 => save=4.063910634805211 => cost: 494.31011527431946exchanging 114-105,10-82 with 114-10,105-82 => save=2.900319415888454 => cost: 491.409795858431Starting path 3/16exchanging 4-13,123-0 with 4-123,13-0 => save=19.13464498164811 => cost: 480.43951154085676exchanging 65-14,52-33 with 65-52,14-33 => save=1.587339787656024 => cost: 478.8521717532007Starting path 4/16exchanging 19-16,55-0 with 19-55,16-0 => save=0.7537655368855951 => cost: 489.3027111852691...exchanging 46-130,39-75 with 46-39,130-75 => save=1.7704562419552907 => cost: 455.4046302803684Starting path 5/16exchanging 37-112,94-0 with 37-94,112-0 => save=11.94699370814547 => cost: 485.7266708931122exchanging 58-110,28-107 with 58-28,110-107 => save=9.304516959883209 => cost: 476.422153933229Starting path 6/16exchanging 62-142,144-0 with 62-144,142-0 => save=23.25445433373693 => cost: 467.79298671639344Starting path 7/16Starting path 8/16exchanging 44-100,48-0 with 44-48,100-0 => save=1.0746379367560621 => cost: 486.7997677569661...exchanging 44-98,8-48 with 44-8,98-48 => save=6.17952602155286 => cost: 476.7280985588525Starting path 9/16Starting path 10/16exchanging 113-5,72-0 with 113-72,5-0 => save=1.4836706701319144 => cost: 484.8140111574118Starting path 11/16Starting path 12/16Starting path 13/16Starting path 14/16Starting path 15/16Starting path 16/16CPU times: user 2.16 s, sys: 40 ms, total: 2.2 sWall time: 2.24 s

Plot improved paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Nice! Several clean paths, but some vans had to do only one distant delivery.

As we see, it is better to increase the max cost again.

%timeresult_paths=sequential_nn(800)

van 1: 792.9220639985251van 2: 793.9102149907774van 3: 749.3577273914822van 4: 769.1224381882283CPU times: user 1.03 s, sys: 18.5 ms, total: 1.05 sWall time: 1.05 s

Plot result paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve them with a 2-opt heuristic

%timeimprove_result_paths()

Starting path 1/4exchanging 84-51,77-0 with 84-77,51-0 => save=11.952005586800027 => cost: 780.9700584117251...exchanging 66-51,147-126 with 66-147,51-126 => save=6.004072981449653 => cost: 732.2204189744828Starting path 2/4exchanging 35-69,29-0 with 35-29,69-0 => save=15.964056951390774 => cost: 777.9461580393865...exchanging 35-129,90-42 with 35-90,129-42 => save=14.532045586175727 => cost: 750.7675349695343Starting path 3/4exchanging 4-124,123-0 with 4-123,124-0 => save=19.23021437016152 => cost: 730.1275130213206...exchanging 65-14,52-33 with 65-52,14-33 => save=1.587339787656024 => cost: 701.2274950765787Starting path 4/4exchanging 19-21,54-0 with 19-54,21-0 => save=12.646642051225456 => cost: 756.4757961370028CPU times: user 3.71 s, sys: 44.8 ms, total: 3.75 sWall time: 3.81 s

Plot improved paths

forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],linewidth=2)plt.scatter(data.x,data.y,c='r',s=30)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Nice, clean paths. But with several overlaps.

Scenario: farthest depot, several vans, KMeans clustering

Compute delivery clusters

fromsklearn.clusterimportKMeanswork_per_van=20#how many deliveries each van will start, ie granularity; a better approach would be capacity.model=KMeans(n_clusters=1+size//work_per_van,init='k-means++',random_state=0)%timemodel.fit(data[['x','y']])

CPU times: user 91.4 ms, sys: 10.6 ms, total: 102 msWall time: 111 msKMeans(algorithm='auto', copy_x=True, init='k-means++', max_iter=300,    n_clusters=8, n_init=10, n_jobs=1, precompute_distances='auto',    random_state=0, tol=0.0001, verbose=0)

Plot clusters

plt.scatter(data.x,data.y,marker='o',c=model.labels_,s=200)plt.scatter(model.cluster_centers_[:,0],model.cluster_centers_[:,1],marker='+',c='r',s=500)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's apply a greedy NN algorithm sequentially, on each of the clusters

%timeresult_paths=greedy_nn_clusters()

van 1: 390.25298877225174van 2: 687.7516428319581van 3: 382.7248611037869van 4: 505.82568809746317van 5: 478.6115152097322van 6: 644.9250908928291van 7: 595.3986361541317van 8: 665.6387979350595CPU times: user 154 ms, sys: 16.6 ms, total: 170 msWall time: 165 ms

Plot result paths

plt.scatter(data.x,data.y,marker='o',c=model.labels_,s=200)plt.scatter(model.cluster_centers_[:,0],model.cluster_centers_[:,1],marker='+',c='r',s=500)forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],'k-',linewidth=2)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Let's improve them with a 2-opt heuristic

%timeimprove_result_paths()

Starting path 1/8exchanging 77-4,83-29 with 77-83,4-29 => save=0.8194748732687884 => cost: 389.43351389898294...exchanging 6-143,128-89 with 6-128,143-89 => save=5.407506824441722 => cost: 324.93008877322006Starting path 2/8exchanging 30-5,34-0 with 30-34,5-0 => save=4.481243727960788 => cost: 683.2703991039973...exchanging 115-95,121-18 with 115-121,95-18 => save=8.322195150021315 => cost: 594.5702251203929Starting path 3/8exchanging 32-59,99-41 with 32-99,59-41 => save=4.114432829335271 => cost: 378.61042827445164exchanging 17-78,26-87 with 17-26,78-87 => save=0.17155006280330554 => cost: 378.43887821164833Starting path 4/8exchanging 42-19,129-106 with 42-129,19-106 => save=3.8728268587745447 => cost: 501.9528612386886...exchanging 37-94,54-0 with 37-54,94-0 => save=1.3747236340181104 => cost: 495.2855372216575Starting path 5/8exchanging 0-147,123-65 with 0-123,147-65 => save=11.990854640529989 => cost: 466.6206605692022...exchanging 135-65,123-0 with 135-123,65-0 => save=3.919229108305899 => cost: 435.93197839715947Starting path 6/8exchanging 56-21,15-0 with 56-15,21-0 => save=6.923848236857282 => cost: 638.0012426559717...exchanging 102-43,118-40 with 102-118,43-40 => save=15.652843812142969 => cost: 602.6373464143627Starting path 7/8exchanging 63-24,107-110 with 63-107,24-110 => save=11.657922849361384 => cost: 583.7407133047703...exchanging 112-36,91-125 with 112-91,36-125 => save=4.722248082275083 => cost: 559.3163570654888Starting path 8/8exchanging 48-111,98-137 with 48-98,111-137 => save=11.5162692633929 => cost: 654.1225286716666...exchanging 8-132,137-98 with 8-137,132-98 => save=7.784370993699362 => cost: 615.4951078634247CPU times: user 2.17 s, sys: 49.8 ms, total: 2.22 sWall time: 2.39 s

Plot improved paths

plt.scatter(data.x,data.y,marker='o',c=model.labels_,s=200)plt.scatter(model.cluster_centers_[:,0],model.cluster_centers_[:,1],marker='+',c='r',s=500)forpathinresult_paths:coords=get_coords(path)plt.plot(coords[0],coords[1],'k-',linewidth=2)plt.scatter(data.iloc[0].x,data.iloc[0].y,marker='*',c='r',s=400)plt.axis([-100,100,-100,100])plt.grid(True)

Very cool subsections! Aproximately the same number of deliveries in each route!

That's it!

Thank you for getting this far, hope you liked! :)