Dijkstra、Bellman-Ford、SPFA、ASP、Floyd-Warshall 演算法分析
圖論中,用來求最短路的方法有很多,適用範圍和時間複雜度也各不相同。
本文主要介紹的演算法的程式碼主要來源如下:
- Dijkstra: Algorithms(《演算法概論》)Sanjoy Dasgupta, Christos Papadimitriou, Umesh Vazirani;《演算法競賽入門經典—訓練指南》劉汝佳、陳峰。
- SPFA (Shortest Path Faster Algorithm): Data Structure and Algorithms Analysis in C, 2nd ed.(《資料結構與演算法分析》)Mark Allen Weiss.
- Bellman-Ford: Algorithms(《演算法概論》)Sanjoy Dasgupta, Christos Papadimitriou, Umesh Vazirani.
- ASP (Acyclic Shortest Paths): Introduction to Algorithms - A Creative Approach(《演算法引論—一種創造性方法》)Udi Manber.
- Floyd-Warshall: Data Structure and Algorithms Analysis in C, 2nd ed.(《資料結構與演算法分析》)Mark Allen Weiss.
它們的使用限制和執行時間如下:
Dijkstra: 不含負權。執行時間依賴於優先佇列的實現,如 O((∣V∣+∣E∣)log∣V∣)
SPFA: 無限制。執行時間O(k⋅∣E∣) (k≪∣V∣)
Bellman-Ford:無限制。執行時間O(∣V∣⋅∣E∣)
ASP: 無圈。執行時間O(∣V∣+∣E∣)
Floyd-Warshall: 無限制。執行時間O(∣V∣^3)
其中 1~4 均為單源最短路徑 (Single Source Shortest Paths) 演算法;5 為全源最短路徑 (All Pairs Shortest Paths) 演算法。順便說一句,為什麼沒有點對點的最短路徑?如果我們只需要一個起點和一個終點,不是比計算一個起點任意終點更節省時間麼?答案還真不是,目前還沒有發現比算從源點到所有點更快的演算法。
分析完時間複雜度之後,再來看平時比賽或者其他使用上的具體情況:
Dijkstra:適用於權值為非負的圖的單源最短路徑,用斐波那契堆的複雜度O(E+VlgV)
BellmanFord:適用於權值有負值的圖的單源最短路徑,並且能夠檢測負圈,複雜度O(VE)
SPFA:適用於權值有負值,且沒有負圈的圖的單源最短路徑,論文中的複雜度O(kE),k為每個節點進入Queue的次數,且k一般<=2,但此處的複雜度證明是有問題的,其實SPFA的最壞情況應該是O(VE).
Floyd:每對節點之間的最短路徑。
這裡給出結論:
(1)當權值為非負時,用Dijkstra。
(2)當權值有負值,且沒有負圈,則用SPFA,SPFA能檢測負圈,但是不能輸出負圈。
(3)當權值有負值,而且可能存在負圈,則用BellmanFord,能夠檢測並輸出負圈。
(4)SPFA檢測負環:當存在一個點入隊大於等於V次,則有負環,後面有證明。
1)圖的表示
本文中,前四個演算法的圖都採用鄰接表表示法,如下:
struct Edge
{
int from;
int to;
int weight;
Edge(int f, int t, int w):from(f), to(t), weight(w) {}
};
int num_nodes;
int num_edges;
vector<Edge> edges;
vector<int> G[max_nodes]; // 每個節點出發的邊編號
int p[max_nodes]; // 當前節點單源最短路中的上一條邊
int d[max_nodes]; // 單源最短路徑長
2)Dijkstra 方法
Dijkstra 方法依據其優先佇列的實現不同,可以寫成幾種時間複雜度不同的演算法。它是圖論-最短路中最經典、常見的演算法。關於這個方法,網上有許多分析,但是我最喜歡的還是《演算法概論》中的講解。為了理解 Dijkstra 方法,首先回顧一下無權最短路的演算法。無權最短路演算法基於 BFS,每次從源點向外擴充套件一層,並且給擴充套件到的頂點標明距離,這個距離就是最短路的長。我們完全可以仿照這個思路,把帶權圖最短路問題規約到無權圖最短路問題——只要把長度大於 1 的邊填充進一些「虛頂點」即可。如下圖所示。
這個辦法雖然可行,但是顯然效率很低。不過,Dijkstra 方法EC,EB,ED分別出發,經過一系列「虛節點」,依次到達D,B,C 。為了不在虛節點處浪費時間,出發之前,我們設定三個鬧鐘,時間分別為4,3,2提醒我們預計在這些時刻會有重要的事情發生(經過實際節點)。更一般地說,假設現在我們處理到了某個頂點u,和u相鄰接的頂點為v1,v2,…,vn,它們和u的距離為d1,d2,…,dn。我們為v1,v2,…,vn各設定一個鬧鐘。如果還沒有設定鬧鐘,那麼設定為d ;如果設定的時間比d晚,那麼重新設定為d(此時我們沿著這條路比之前的某一條路會提前趕到)。每次鬧鐘響起,都說明可能經過了實際節點,我們都會更新這些資訊,直到不存在任何鬧鐘。綜上所述,也就是隨著 BFS 的進行,我們一旦發現更近的路徑,就立即更新路徑長,直到處理完最後(最遠)的一個頂點。由此可見,由於上述「虛頂點」並非我們關心的實際頂點,因此 Dijkstra 方法的處理方式為:直接跳過了它們。
還需要解決的一個問題,就是鬧鐘的管理。鬧鐘一定是從早到晚按順序響起的,然而我們設鬧鐘的順序卻不一定按照時間升序,因此需要一個優先佇列來管理。Dijkstra 方法實現的效率嚴重依賴於優先佇列的實現。一個使用標準庫容器介面卡 priority_queue
的演算法版本如下:
typedef pair<int, int> HeapNode;
void Dijkstra(int s)
{
priority_queue< HeapNode, vector<HeapNode>, greater<HeapNode> > Q;
for (int i=0; i<num_nodes; ++i)
d[i] = __inf;
d[s] = 0;
Q.push(make_pair(0, s));
while (!Q.empty()) {
pair<int, int> N = Q.top();
Q.pop();
int u = N.second;
if (N.first != d[u]) continue;
for (int i=0; i<G[u].size(); ++i) {
Edge &e = edges[G[u][i]];
if (d[e.to] > d[u] + e.weight) {
d[e.to] = d[u] + e.weight;
p[e.to] = G[u][i];
Q.push(make_pair(d[e.to], e.to));
}
}
}
}
3)Bellman-Ford:一個簡單的想法
Dijkstra 方法的本質是進行一系列如下的更新操作:
d(v)=min{d(v), d(u)+l(u,v)}
然而,如果邊權含有負值,那麼 Dijkstra 方法將不再適用。原因解釋如下。
假設最終的最短路徑為:
s→u1→u2→u3→…→uk→t
不難看出,如果按照 (s, u1), (u1, u2), …,(uk, t) 的順序執行上述更新操作,最終的最短路徑一定是正確的。而且,只要保證上述更新操作全部按順序執行即可,並不要求上述更新操作是連續進行的。Dijkstra 演算法所執行的更新序列是經過選擇的,而選擇基於這一假設:s→t的最短路一定不會經過和距離大於l(s, t)的點。對於正權圖這一假設是顯然的,對於負權圖這一假設是錯誤的。
因此,為了求出負權圖的最短路徑,我們需要保證一個合理的更新序列。但是,我們並不知道最終的最短路徑!因此一個簡單的想法就是:更新所有的邊,每條邊都更新∣V∣−1次。由於多餘的更新操作總是無害的,因此演算法(幾乎)可以正確執行。等等,為什麼是∣V∣−1次?這是由於,任何含有∣V∣個頂點的圖兩個點之間的最短路徑最多含有∣V∣−1條邊。這意味著最短路不會包含環。理由是,如果是負環,最短路不存在;如果是正環,去掉後變短;如果是零環,去掉後不變。
演算法實現中唯一一個需要注意的問題就是負值圈 (negative-cost cycle)。負值圈指的是,權值總和為負的圈。如果存在這種圈,我們可以在裡面滯留任意長而不斷減小最短路徑長,因此這種情況下最短路徑可能是不存在的,可能使程式陷入無限迴圈。好在,本文介紹的幾種演算法都可以判斷負值圈是否存在。對於 Bellman-Ford 演算法來說,判斷負值圈存在的方法是:在次迴圈之後再執行一次迴圈,如果還有更新操作發生,則說明存在負值圈。
Bellman-Ford 演算法的程式碼如下:
bool Bellman_Ford(int s)
{
for (int i=0; i<num_nodes; ++i)
d[i] = __inf;
d[s] = 0;
for (int i=0; i<num_nodes; ++i) {
bool changed = false;
for (int e=0; e<num_edges; ++e) {
if (d[edges[e].to] > d[edges[e].from] + edges[e].weight
&& d[edges[e].from] != __inf) {
d[edges[e].to] = d[edges[e].from] + edges[e].weight;
p[edges[e].to] = e;
changed = true;
}
}
if (!changed) return true;
if (i == num_nodes && changed) return false;
}
return false; // 程式應該永遠不會執行到這裡
}
註記:
- 如果某次迴圈沒有更新操作發生,以後也不會有了。我們可以就此結束程式,避免無效的計算。
- 上述程式中第 11 行的判斷:如果去掉這個判斷,只要圖中存在負值圈函式就會返回
false
。否則,僅在給定源點可以達到負值圈時才返回false
。
4)SPFA:一個改進的想法
看來,Bellman-Ford 演算法多少有些浪費。這裡介紹的 SPFA 可以算作是 Bellman-Ford 的佇列改進版。在 Dijkstra 方法中,隨著 BFS 的進行,最短路一點一點地「生長」。然而如果存在負權,我們的演算法必須允許「繞回去」的情況發生。換句話說,我們需要在某些時候撤銷已經形成的最短路。同時,我們還要改變 Bellman-Ford 演算法盲目更新的特點,只更新有用的節點。SPFA 中,一開始,我們把源點 s放入佇列中,然後每次迴圈讓一個頂點u出隊,找出所有與u鄰接的頂點v,更新最短路,並當v不在佇列裡時將它入隊。迴圈直到佇列為空(沒有需要更新的頂點)。
可以顯示出 SPFA 和 Bellman-Ford 演算法相比的一個重大改進的最經典的一個例子,就是圖為一條鏈的情形。
容易看出,如果存在負值圈,這個演算法將無限迴圈下去。因此需要一個機制來確保演算法得以中止。由於最短路最長只含有∣V∣−1條邊,因此如果任何一個頂點已經出隊∣V∣+1次,演算法停止執行。
證明如果有負環當且僅當存在一個點入佇列次數大於等於V次:
對於某個點v,我們已知s到v的鬆弛路徑的邊的數量最多為V-1。
我這裡說的鬆弛路徑指的是:比如s直接鬆弛v,這樣就有一條鬆弛路徑:s->v 。s鬆弛a,a鬆弛v,則s->a->v就是一條鬆弛路徑。
對於所有s到v的鬆弛路徑來說,當鬆弛路徑邊的數量相等時,v只入隊一次。
比如有鬆弛路徑:
s->a->x->v
s->b->x->v
s->c->z->v,可以看出v只入隊一次。
因為s到v的鬆弛路徑的長度最多可以有V-1種變化,所以v最多入隊V-1次。
舉個例子:
假設有一個圖,點集為{s,a,b,c,v},則最多可能的鬆弛路徑有:
s->v
s->a->v
s->b->v
s->c->v
s->a->b->v
s->a->c->v
s->b->c->v
s->a->b->c->v
則鬆弛路徑的邊數變化有1,2,3,4,所以v入隊為4次,即V-1次。
所以我們可以說每個點最多入隊V-1次,因此我們求最壞情況為每個點都入隊V-1次,所以此時:
這裡舉個最壞情況的例子。
當然我們可能考慮,當給定一個V的值,E的值,比如E=2V,怎麼給出一個圖,使得對此圖執行SPFA演算法的複雜度為O(VE).
我們這裡假定圖是連通的,所以E>=V-1。
方法如下:
(1)我們首先將圖組成一個鏈,即如下圖所示:
這樣就用去了V-1條邊。
(2)分別新增v0連向v2,v3,….vk的邊,我們要新增的這些邊的權值要滿足v0先更新vk,v0更新vk-1後vk-1還能更新vk,以此類推,如下圖所示:
(3)以vk,vk-1,…..v1的順序新增權值為正無窮的自環,且不斷迴圈,這個步驟是為了保持v1到vk點的出度保持一致,所以這樣做,如下圖:
這樣我們就構造了一個能夠讓SPFA跑出O(VE)的圖了,原因如下:
因為我們E的值和V的值是不確定的,所以很有可能不能夠完成上述的這些構造,我們會分析當沒有剩餘的邊構造上面的步驟(2)時的複雜度(也就是說E<=2V-3,因為第一步連成一個鏈需要V-1條邊,而第二步v0連出去的邊需要V-2),和有足夠的邊能構造上面的圖這兩種情況。
(1)如果E<=2V-3
因為E<=2V-3,所以E=O(V),所以只要能夠求出複雜度是O(V^2),即可說為O(VE).
我們要計算所有點的入隊次數和訪問的邊數。
V0出度為 E-V+2,V0入隊1次。
V1出度為1,V1入隊為1次。
V2出度為1,V2入隊為2次(分別為v0鬆弛v2,v1鬆弛v2)。
V3出度為1,V3入隊為3次。
….
V(e-v+2)出度為1,入隊次數為(E-V+2)次。
後面V(e-v+3),V(e-v+4),…..V(k-1)的出度為1,入隊次數為E-V+2次。 這些點的個數為V-(E-V+3)-1 = 2V-E-4。
vk出度為0,vk入隊為E-V+2次。
所以總共的訪問的邊數為:
(2)如果E>2V-3
此時構造圖的第二步已經完畢,所以後面剩餘的邊只需要不斷新增自環保持出度平衡即可。
V0出度為V-1,入隊1次。
V1到Vk出度為(E-V+2)/(V-1) 或(E-V+2)/(V-1)+1。
v1到vk的入隊次數分別為1,2,3,…..V-1。
所以總共訪問邊數為:
SPFA 的程式碼如下:
bool in_queue[max_nodes];
int cnt[max_nodes];
bool SPFA(int s)
{
int u;
queue<int> Q;
memset(in_queue, 0, sizeof(in_queue));
memset(cnt, 0, sizeof(cnt));
for (int i=0; i<num_nodes; ++i)
d[i] = __inf;
d[s] = 0;
in_queue[s] = true;
Q.push(s);
while (!Q.empty()) {
in_queue[u=Q.front()] = false;
Q.pop();
for (int i=0; i<G[u].size(); ++i) {
Edge &e = edges[G[u][i]];
if (d[e.to] > d[u] + e.weight) {
d[e.to] = d[u] + e.weight;
p[e.to] = G[u][i];
if (!in_queue[e.to]) {
Q.push(e.to);
in_queue[e.to] = true;
if (++cnt[e.to] > num_nodes)
return false;
}
}
}
}
return true;
}
我們已經給出 SPFA 的執行時間為O(k⋅∣E∣) (k≪∣V∣)。實際上,可以期望k<2。但是可以構造出使 SPFA 演算法變得很慢的針對性資料。
5)Acyclic Shortest Path
如果圖是無圈的(acyclic)(或稱為有向無環圖,DAG),那麼情況就變的容易了。我們可以構造一個拓撲升序序列,由拓撲排序的性質,無圈圖的任意路徑中,頂點都是沿著「拓撲升序序列」排列的,因此我們只需要按照這個序列執行更新操作即可。顯然,這樣可以線上性時間內解決問題。
實現上,拓撲排序和更新可以一趟完成。這種演算法的程式碼如下:
int indegree[max_nodes];
void asp(int s)
{
queue<int> Q;
for (int i=0; i<num_nodes; ++i) {
d[i] = __inf;
indegree[i] = 0;
}
for (int i=0; i<num_edges; ++i)
++indegree[edges[i].to];
for (int i=0; i<num_nodes; ++i)
if (indegree[i] == 0) Q.push(i);
d[s] = 0;
while (!Q.empty()) {
int w = Q.front();
Q.pop();
for (int i=0; i<G[w].size(); ++i) {
if (d[edges[G[w][i]].to] > d[w] + edges[G[w][i]].weight && d[w] != __inf) {
d[edges[G[w][i]].to] = d[w] + edges[G[w][i]].weight;
p[edges[G[w][i]].to] = G[w][i];
}
if (–indegree[edges[G[w][i]].to] == 0)
Q.push(edges[G[w][i]].to);
}
}
}
6)Floyd-Warshall
與前面四種不同,Floyd-Warshall 演算法是所謂的「全源最短路徑」,也就是任意兩點間的最短路徑。它並不是對單源最短路徑
∣V∣次迭代的一種漸進改進,但是對非常稠密的圖卻可能更快,因為它的迴圈更加緊湊。而且,這個演算法支援負的權值。
Floyd-Warshall 演算法如下:
int dist[max_nodes][max_nodes]; // 記錄路徑長
int path[max_nodes][max_nodes]; // 記錄實際路徑
bool Floyd_Warshall()
{
for (int i=0; i<num_nodes; ++i)
for (int j=0; j<num_nodes; ++j)
path[i][j] = j;
for (int k=0; k<num_nodes; ++k) {
for (int i=0; i<num_nodes; ++i) {
for (int j=0; j<num_nodes; ++j) {
if (dist[i][j] > dist[i][k] + dist[k][j]
&& dist[i][k] != __inf && dist[k][j] != __inf) {
dist[i][j] = dist[i][k] + dist[k][j];
path[i][j] = path[i][k];
if (i == j && dist[i][j] < 0)
return false;
}
}
}
}
return true;
}
其中 dist
陣列應初始化為鄰接矩陣。需要提醒的是, dist[i][i]
實際上表示「從頂點 i 繞一圈再回來的最短路徑」,因此圖存在負環當且僅當 dist[i][i]<0。初始化時, dist[i][i]
可以初始化為0
,也可以初始化為∞。
7)顯示實際路徑
顯示實際路徑實際上非常簡單。利用前四個演算法提供的 int *p
,還原實際路徑的一個方法如下:
void printpath(int from, int to, bool firstcall = true)
{
if (from == to) {
printf(”%d”, from);
return;
}
if (p[to] == -1) return;
if (firstcall) printf(“%d ->”, from);
int v = edges[p[to]].from;
if (v == from) {
if (firstcall) printf(“ %d”, to);
return;
}
printpath(from, v, false);
printf(” %d ->”, v+1);
if (firstcall) printf(“ %d”, to);
利用 Floyd-Warshall 演算法提供的 int **path
,還原實際路徑的一個方法如下:
void showpath(int from, int to)
{
int c = from;
printf(”%d -> %d:(%2d) %d”, from, to, dist[from][to], from);
while (c != to) {
printf(” -> %d”, path[c][to]);
c = path[c][to];
}
printf(”\n”);
}