基于Gurobi軟件Callback功能的旅行商問題求解

度巍 陳昊澤

摘要：作為經典組合優化問題，旅行商問題（Traveling Salesman Problem簡稱TSP） 一直是大學交通運輸與應用數學等專業的教學與科研熱點。在基于混合整數規劃模型的TSP求解中，需要解決如何避免出現子環路問題，Gurobi作為當前最先進的運籌優化軟件，其具有的Callback功能使模型在求解過程中，動態地添加子環路約束成為可能。文章針對當前相關網絡資源存在的問題，構建了用Python編寫的基于Callback功能動態添加子環路消除約束的TSP求解代碼，通過多個算例驗證了代碼的求解可行性，為逐步將Gurobi引入課堂教學提供了素材。

關鍵詞：旅行商問題;子環路消除;Gurobi;Callback功能

中圖分類號：G642? ? ? ? 文獻標識碼：A

文章編號：1009-3044（2022）25-0009-02

開放科學（資源服務） 標識碼（OSID） ：

1 引言

新一代大規模優化軟件Gurobi，因其在求解數學規劃問題的卓越性能，逐漸成為高校師生運籌學教學與科研首選軟件。在Gurobi的高級功能中，Callback函數可獲取求解過程中的信息，動態加入新的約束條件或者其他算法，為實現各種復雜問題的求解創造了條件。本文通過Callback功能，在求解旅行商問題過程中，動態添加子環路消除約束條件，實現了一種新的求解旅行商問題方法，并通過與現有模型在不同問題規模求解時間的對比，幫助學生加深對旅行商問題的理解，為開展Gurobi軟件的教學與相關科研提供了素材。

2 旅行商問題的數學規劃模型

組合優化領域中的旅行商問題（Traveling Salesman Problem簡稱TSP） ，廣泛應用物流配送、電纜和光纜布線等領域[1]，如何找到大規模TSP的最優解一直是國內外學者的研究熱點[2]，近幾十年來涌現出各種求解算法。旅行商問題一般表述為：某旅行推銷員要去若干個城市推銷商品，然后回到其出發城市，已知任意兩個城市之間的行走距離，求出該推銷員經過每個城市有且一次同時總距離最短的巡回路線。若城市個數為[n]，城市集合為[V={1，2，…，n}]，從城市[i]到城市[j]的距離為[cij]，若在巡回路線中，從城市[i]直接走到城市[j]，0-1變量[xij]取值為1，否則取值為0，則求解旅行商問題涉及如下整數規劃模型：[3]

[min Z=i=1nj=1ncijxij? ?i，j∈V，i≠js.t.? ? ? ? ? ?i=1nxij=1，? ??j∈V，i≠j? ? ? ? ? ? ? ? j=1nxij=1，? ??i∈V，i≠j? ? ? ? ? ? ? ? xij∈{0，1}， ?i，j∈V，i≠j]? ? ? ? （1）

模型（1）中的目標函數表示巡回路線總的長度，第1類與第2類約束條件分別表示每個城市節點只能進出一次。然而如果僅僅是求解模型（1），獲得的解有可能出現子環路情形，即解所對應的路線中存在少于[n]的若干城市節點首尾相連的情形。為了避免該子環路問題，有學者做過大量的研究，如在模型（1）上增加如下MTZ約束：

[ui-uj+nxij≤n-1， ?i，j∈V，i≠j]? ? ? ? （2）

（2）式中，實數決策變量[ui]表示城市[i]在巡回路線中的序號。添加MTZ約束是目前大多數運籌優化軟件在避免出現子環路時采取的方法。如高校運籌學教學應用廣泛的LINGO軟件，在其模型庫中，求解TSP采取的是文獻[4]對應的增強MTZ條件：

[uj≥ui+xij-（n-2）（1-xij）+（n-3）xji? ?i，j∈{2…n}，i≠jui≤n-1-（n-2）x1i? ? ? ? ? ? ? ? ? ? ? ? ? ? ??i∈{2…n}ui≥1+（n-2）x1i? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ??i∈{2…n}]? ? ? （3）

解決避免出現子環路的另一個思路是在模型（1）上添加如下子環路消除約束：

[i，j∈Sxij≤S-1， 2≤S≤n-1，S?V]? ? ? （4）

其中[S]表示所有的子環路集合，[S]表示[S]包含的城市節點個數，（4）式的理解非常直觀，即讓[S]中依次相連的節點之間對應的[xij]之和少于[S]，從而把子環路出現時對應的解排除掉。然而，由于該約束數量與[2n]同階，過去很難應用于實際問題的求解中。目前Gurobi軟件具有的Callback功能[5]，可以動態地將求解過程中出現的每個子環路所對應的消除約束（4）式依次添加到模型中，基于Callback功能的完整求解流程如圖1所示：

求解TSP流程圖

3 模型求解的關鍵代碼

基于上述途徑求解TSP的代碼用Python語言實現，通過導入gurobipy庫實現調用Gurobi的Callback功能，部分代碼參考了網絡資源[6]。筆者在對文獻[6]的考察中，發現其相關代碼并不能求出最優解，進一步對代碼做了改進完善。由于代碼較長，只給出關鍵代碼部分。

檢查當前解所對應的巡游路線函數段：

def subtour（graph）：

unvisited = range（0，n）? #構造探尋可能構成子環路城市節點序列

cycle = range（0， n）? # cycle用來返回找到的子環路

edges = findEdges（graph）

edges = tuplelist（edges）

tempcycle=[]

thiscycle1 = [unvisited[0]] #從unvisited的第一個節點可以尋找子環路

unvisited.remove（0）

while len（thiscycle1）

point1=thiscycle1[-1]

neighbors1 = [j for i， j in edges.select（point1， '*'） if j in unvisited]

if? len（neighbors1） >= 1：

thiscycle1.extend（neighbors1）

if （thiscycle1[-1]， thiscycle1[0]） in edges：

tempcycle=thiscycle1

break

cycle = tempcycle

return cycle

對當前解是否存在子環路進行檢查，如果存在子環路，構造出對應的約束條件（4）函數段：

def subtourelim（model， where）：

if （where == GRB.Callback.MIPSOL）：

x_value = np.zeros（[nodeNum ， nodeNum ]）

for m in model.getVars（）：

if （m.varName.startswith（'x'））：

a = （int）（m.varName.split（'_'）[1]）

b = （int）（m.varName.split（'_'）[2]）

x_value[a][b] = model.cbGetSolution（m）

tour = subtour（x_value）? #從subtour函數獲取當前解對應的環路

if （len（tour） < n）：? #如果得到的環路包含節點個數小于城市總數，即為子環路

# 將該子環路所對應的子環路消除約束添加到模型中

tt1 = LinExpr（0）

for i in range（0，len（tour））：

if （i < len（tour）-1）：

tt1.addTerms（1， X[tour[i]，tour[i+1]]）

elif （i == len（tour）-1）：

tt1.addTerms（1， X[tour[i]，tour[0]]）

model.cbLazy（tt1<= len（tour） - 1）

不同于文獻[6]在構建基本模型時額外添加一個虛擬起始點的做法，本文代碼直接構建模型（1）中的決策變量與目標函數：

model = Model（'TSP'）

X = {}? #構建決策變量

mu = {}

for j in range（nodeNum）：

if （i ！= j）：

X[i， j] = model.addVar（vtype=GRB.BINARY， name='x_' + str（i） + '_' + str（j））

obj = LinExpr（0）? # 構建目標函數

for key in X.keys（）：

i = key[0]

j = key[1]

if （i < nodeNum and j < nodeNum）： obj.addTerms（cost[key[0]][key[1]]， X[key]）

4 實例求解分析

結合上面三個關鍵段以及導入數據，調用Gurobi求解模型等代碼部分，能實現對TSP的求解。在Windows10操作系統、8G內存、Inter Core2.5GHz實驗平臺上，做了一系列不同規模TSP的求解實驗，對較小規模的TSP，本文構造的模型可以較快地求出最優巡回路徑，如針對文獻[6]提到的Solomon算例集中的C201節點數據，取前11個節點，能在幾秒內計算出最優巡游序列為1→7→9→10→11→6→3→2→8→4→5→1，進一步將本文的模型與基于MTZ約束的模型在不同節點數量下的相同問題求解時間進行了對比，結果如表1所示：

可以看到，在實例規模較小時，兩種模型求解時間差別不大，但隨著問題規模的增加，本文構建的Callback動態添加子環路消除約束模型在求解時間上迅速增加。這是因為子環路個數隨著問題規模的增加呈現指數增加，模型在每一次添加新的子環路消除約束后又要重新計算新的解，而相比較下，基于MTZ約束的模型一次性添加好所有約束，在求解時間上明顯更有優勢，這很好解釋了當前類似LINGO、GAMS等軟件在求解TSP時基于MTZ約束的原因。

5 結束語

旅行商問題是高校交通運輸本科專業開設的諸如《交通運籌學》《交通系統分析》等課程中重要知識點，目前課堂的教學往往是對相關模型與求解算法進行介紹，這很難使學生對該問題形成深刻印象，更不能體會到求解旅行商問題的難度與應用價值。本文基于Gurobi軟件的求解優化問題強大功能，實現了過去無法做到的通過添加子環路消除約束求解旅行商問題途徑，讓學生直觀感受到不同的求解模型在計算時間上的差別，加深了對該問題的理解，直觀體會到問題的計算復雜程度。同時也將本科生在計算機課程學到的Python語言應用于專業知識實踐，實現了大學各個課程間的相互融入，為學生進一步學習車輛路徑問題、交通分配模型打下基礎。

參考文獻：

[1] William J.Cook.迷茫的旅行商一個無處不在的計算機算法問題[M]. 隋春寧，譯.北京：人民郵電出版社，2013.

[2] 李汝佳.基于蟻群算法的旅游路線規劃問題研究[J].電腦知識與技術，2019，15（8）：137-140.

[3] 韓中庚.運籌學及其工程應用[M].北京：清華大學出版社，2014.

[4] Desrochers M，Laporte G.Improvements and extensions to the Miller-Tucker-Zemlin subtour elimination constraints[J].Operations Research Letters，1991，10（1）：27-36.

[5] Gurobi OptimizationCompany. Gurobi documents[EB/OL].[2021-10-01] https：//www.gurobi.com/documentation/9.5/refman/index.html

[6] 劉興祿.TSP中兩種不同消除子環路的方法及callback實現[EB/OL].[2021-01-16]. https：//blog.csdn.net/weixin_398332 90/article/details/113452735？utm_medium=distribute.pc_relevant.none-task-blog-2～default～baidujs_title～default-0.no_ search_link&spm=1001.2101.3001.4242.1.

【通聯編輯：王力】