求單源最短路的SPFA算法的全稱是:Shortest Path Faster Algorithm。
SPFA算法是西南交通大學段凡丁于1994年發表的。
從名字我們就可以看出,這種算法在效率上一定有過人之處。
很多時候,給定的圖存在負權邊,這時類似Dijkstra等算法便沒有了用武之地,而Bellman-Ford算法的複雜度又過高,SPFA算法便派上用場了。有人稱spfa算法是最短路的萬能算法。
簡潔起見,我們約定有向權重圖G不存在負權回路,即最短路徑一定存在。當然,我們可以在執行該算法前做一次拓撲排序,以判斷是否存在負權回路。
我們用數組dis記錄每個結點的最短路徑估計值,可以用鄰接矩陣或鄰接表來存儲圖G,推薦使用鄰接表。
spfa的算法思想(動态逼近法):
設立一個先進先出的隊列q用來儲存待優化的結點,優化時每次取出隊首結點u,并且用u點目前的最短路徑估計值對離開u點所指向的結點v進行松弛操作,如果v點的最短路徑估計值有所調整,且v點不在目前的隊列中,就将v點放入隊尾。這樣不斷從隊列中取出結點來進行松弛操作,直至隊列空為止。
松弛操作的原理是著名的定理:“三角形兩邊之和大于第三邊”,在資訊學中我們叫它三角不等式。所謂對結點i,j進行松弛,就是判定是否dis[j]>dis[i]+w[i,j],如果該式成立則将dis[j]減小到dis[i]+w[i,j],否則不動。
下面舉一個執行個體來說明SFFA算法是怎樣進行的:
和廣搜bfs的差別:
SPFA 在形式上和廣度(寬度)優先搜尋非常類似,不同的是bfs中一個點出了隊列就不可能重新進入隊列,但是SPFA中一個點可能在出隊列之後再次被放入隊列,也就是一個點改進過其它的點之後,過了一段時間可能本身被改進(重新入隊),于是再次用來改進其它的點,這樣反複疊代下去。
算法的描述:
void spfa(s); //求單源點s到其它各頂點的最短距離
for i=1 to n do { dis[i]=∞; vis[i]=false; } //初始化每點到s的距離,不在隊列
dis[s]=0; //将dis[源點]設為0
vis[s]=true; //源點s入隊列
head=0; tail=1; q[tail]=s; //源點s入隊, 頭尾指針賦初值
while head<tail do {
head+1; //隊首出隊
v=q[head]; //隊首結點v
vis[v]=false; //釋放對v的标記,可以重新入隊
for 每條邊(v,i) //對于與隊首v相連的每一條邊
if (dis[i]>dis[v]+a[v][i]) //如果不滿足三角形性質
dis[i] = dis[v] + a[v][i] //松弛dis[i]
if (vis[i]=false) {tail+1; q[tail]=i; vis[i]=true;} //不在隊列,則加入隊列
}
最短路徑本身怎麼輸出?
在一個圖中,我們僅僅知道結點A到結點E的最短路徑長度,有時候意義不大。這個圖如果是地圖的模型的話,在算出最短路徑長度後,我們總要說明“怎麼走”才算真正解決了問題。如何在計算過程中記錄下來最短路徑是怎麼走的,并在最後将它輸出呢?
我們定義一個path[]數組,path[i]表示源點s到i的最短路程中,結點i之前的結點的編号(父結點),我們在借助結點u對結點v松弛的同時,标記下path[v]=u,記錄的工作就完成了。
如何輸出呢?我們記錄的是每個點前面的點是什麼,輸出卻要從最前面到後面輸出,這很好辦,遞歸就可以了:
c++ code:
void printpath(int k){
if (path[k]!=0) printpath(path[k]);
cout << k << ' ';
}
pascal code:
procedure printpath(k:longint);
begin
if path[k]<>0 then printpath(path[k]);
write(k,' ');
end;
spfa算法模闆(鄰接矩陣):
c++ code:
void spfa(int s){
for(int i=0; i<=n; i++) dis[i]=99999999; //初始化每點i到s的距離
dis[s]=0; vis[s]=1; q[1]=s; 隊列初始化,s為起點
int i, v, head=0, tail=1;
while (head<tail){ 隊列非空
head++;
v=q[head]; 取隊首元素
vis[v]=0; 釋放隊首結點,因為這節點可能下次用來松弛其它節點,重新入隊
for(i=0; i<=n; i++) 對所有頂點
if (a[v][i]>0 && dis[i]>dis[v]+a[v][i]){
dis[i] = dis[v]+a[v][i]; 修改最短路
if (vis[i]==0){ 如果擴充結點i不在隊列中,入隊
tail++;
q[tail]=i;
vis[i]=1;
}
}
}
}
pascal code:
procedure spfa(s:longint);
var i,j,v,head,tail:longint;
begin
for i:=0 to n do dis[i]:=99999999;
dis[s]:=0; vis[s]:=true; q[1]:=s;
head:=0;tail:= 1;
while head<tail do
begin
inc(head);
v:=q[head];
vis[v]:=false;
for i:=0 to n do
if dis[i]>dis[v]+a[v,i] then
begin
dis[i]:= dis[v]+a[v,i];
if not vis[i] then
begin
inc(tail);
q[tail]:=i;
vis[i]:=true;
end;
end;
end;
end;
【程式1】暢通工程 (laoj1138)
某省自從實行了很多年的暢通工程計劃後,終于修建了很多路。不過路多了也不好,每次要從一個城鎮到另一個城鎮時,都有許多種道路方案可以選擇,而某些方案要比另一些方案行走的距離要短很多。這讓行人很困擾。
現在,已知起點和終點,請你計算出要從起點到終點,最短需要行走多少距離。
輸入格式:
第一行包含兩個正整數N和M(0<N<200,0<M<1000),分别代表現有城鎮的數目和已修建的道路的數目。城鎮分别以0~N-1編号。
接下來是M行道路資訊。每一行有三個整數A,B,X(0<=A,B<N,A!=B,0<X<10000),表示城鎮A和城鎮B之間有一條長度為X的雙向道路。
再接下一行有兩個整數S,T(0<=S,T<N),分别代表起點和終點。
輸出格式:
輸出最短需要行走的距離。如果不存在從S到T的路線,就輸出-1。
樣例輸入1:
3 3
0 1 1
0 2 3
1 2 1
0 2
樣例輸入2:
3 1
0 1 1
1 2
樣例輸出1:
2
樣例輸出2:
-1
【分析】 注意本題可能有結點為0的頂點,我就在這上面wa了很多次。并且有可能兩個城鎮之間有多條道路,我們要保留最小的那條。
pascal code(鄰接矩陣):
var i,n,m,s,t,x,y,z:longint; s:起點;t:終點
a,b:array[0..201,0..201] of longint; b[x,c]存與x相連的第c個邊的另一個結點y
q:array[0..10001] of integer; 隊列
vis:array[0..201] of boolean; 是否入隊的标記
dis:array[0..201] of longint; 到起點的最短路
procedure spfa(s:longint);
var i,j,v,head,tail:longint;
begin
fillchar(q,sizeof(q),0);
fillchar(vis,sizeof(vis),false);
for i:=0 to n do dis[i]:=99999999;
dis[s]:=0; vis[s]:=true; q[1]:=s; 隊列的初始狀态,s為起點
head:=0;tail:= 1;
while head<tail do 隊列不空
begin
inc(head);
v:=q[head]; 取隊首元素
vis[v] := false; 釋放結點,一定要釋放掉,因為這節點有可能下次用來松弛其它節點
for i:=1 to b[v,0] do
if dis[b[v,i]]>dis[v]+a[v,b[v,i]] then
begin
dis[b[v,i]]:=dis[v]+a[v,b[v,i]]; 修改最短路
if not vis[b[v,i]] then 擴充結點入隊
begin
inc(tail);
q[tail]:=b[v,i];
vis[b[v,i]]:=true;
end;
end;
end;
end;
begin
read(n, m); //n結點數;m邊數
fillchar(a,sizeof(a),0);
for i:=1 to m do
begin
readln(x,y,z); x,y一條邊的兩個結點;z這條邊的權值
if (a[x,y]<>0)and(z>a[x,y]) then continue;如果兩頂點間有多條邊,保留最小的一條
inc(b[x,0]);b[x,b[x,0]]:=y;a[x,y]:=z; b[x,0]以x為一個結點的邊的條數
inc(b[y,0]);b[y,b[y,0]]:=x;a[y,x]:=z;
end;
readln(s,t); 讀入起點與終點
spfa(s);
if dis[t]<>99999999 then writeln(dis[t]) else writeln(-1);
end.
C++ code(鄰接矩陣):
#include <iostream>
using namespace std;
int q[10001], dis[201], a[201][201], b[201][201];
bool vis[201];
int n, m, s, t;
void spfa(int s){
for(int i=0; i<=n; i++) dis[i]=99999999;
dis[s]=0; vis[s]=1; q[1]=s; 隊列的初始狀态,s為起點
int i, v, head=0, tail=1;
while (head<tail){ 隊列不空
head++;
v=q[head]; 取隊首元素
vis[v]=0; 釋放結點,一定要釋放掉,因為這節點有可能下次用來松弛其它節點
for(i=1; i<=b[v][0]; i++)
if (dis[b[v][i]] > dis[v]+a[v][b[v][i]]){
dis[b[v][i]] = dis[v]+a[v][b[v][i]]; 修改最短路
if (vis[b[v][i]]==0){ 擴充結點入隊
tail++;
q[tail]=b[v][i];
vis[b[v][i]]=1;
}
}
}
}
int main(){
int x, y, z;
cin >> n >> m; //n結點數;m邊數
for(int i=0; i<m; i++){
cin >> x >> y >> z; x,y一條邊的兩個結點;z這條邊的權值
if (a[x][y]!=0 && z>a[x][y]) continue;如果兩頂點間有多條邊,保留最小的一條
b[x][0]++; b[x][b[x][0]]=y; a[x][y]=z; b[x,0]以x為一個結點的邊的條數
b[y][0]++; b[y][b[y][0]]=x; a[y][x]=z;
}
cin >> s >> t; 讀入起點與終點
spfa(s);
if (dis[t]!=99999999) cout << dis[t] << endl;
else cout << -1 << endl;
return 0;
}
spfa優化——深度優先搜尋dfs
在上面的spfa标準算法中,每次更新(松弛)一個結點u時,如果該結點不在隊列中,那麼直接入隊。
但是有負環時,上述算法的時間複雜度退化為O(nm)。能不能改進呢?
那我們試着使用深搜,核心思想為
每次從更新一個結點u時,從該結點開始遞歸進行下一次疊代。
使用dfs優化spfa算法:
pascal code:
procedure spfa(s:longint);
var i:longint;
begin
for i:=1 to b[s,0] do //b[s,0]是從頂點s發出的邊的條數
if dis[b[s,i]]>dis[s]+a[s,b[s,i]] then //b[s,i]是從s發出的第i條邊的另一個頂點
begin
dis[b[s,i]]:=dis[s]+a[s,b[s,i]];
spfa(b[s,i]);
end;
end;
C++ code:
void spfa(int s){
for(int i=1; i<=b[s][0]; i++) //b[s,0]是從頂點s發出的邊的條數
if (dis[b[s][i]>dis[s]+a[s][b[s][i]]){ //b[s,i]是從s發出的第i條邊的另一個頂點
dis[b[s][i]=dis[s]+a[s][b[s][i]];
spfa(b[s][i]);
}
}
相比隊列,深度優先搜尋有着先天優勢:在環上走一圈,回到已周遊過的結點即有負環。絕大多數情況下的時間複雜度為O(m)級别。
那我們試着使用深搜,核心思想為
每次從更新一個結點u時,從該結點開始遞歸進行下一次疊代。
對于WorldRings(ACM-ICPC Centrual European 2005)這道題,676個點,100000條邊,查找負環dfs僅僅需219ms。
一個簡潔的資料結構和算法在一定程度上解決了大問題。
判斷存在負環的條件:重新經過某個在目前搜尋棧中的結點。
【程式1】暢通工程 laoj1138 spfa算法(dfs):
pascal code:
var i,n,m,s,t,x,y,z:longint;
a,b:array[0..201,0..201] of longint;
q:array[0..10001] of integer;
vis:array[0..201] of boolean;
dis:array[0..201] of longint;
procedure spfa(s:longint);
var i:longint;
begin
for i:=1 to b[s,0] do
if dis[b[s,i]]>dis[s]+a[s,b[s,i]] then
begin
dis[b[s,i]]:=dis[s]+a[s,b[s,i]];
spfa(b[s,i]);
end;
end;
begin
read(n, m);
fillchar(a,sizeof(a),0);
for i:=1 to m do
begin
readln(x,y,z);
if (a[x,y]<>0)and(z>a[x,y]) then continue;
inc(b[x,0]);b[x,b[x,0]]:=y;a[x,y]:=z;
inc(b[y,0]);b[y,b[y,0]]:=x;a[y,x]:=z;
end;
readln(s,t);
for i:=0 to n do dis[i]:=99999999;
dis[s]:=0;
spfa(s);
if dis[t]<>99999999 then writeln(dis[t]) else writeln(-1);
end.
C++ code:
#include <iostream>
using namespace std;
int q[10001], dis[201], a[201][201], b[201][201];
bool vis[201];
int n, m, s, t;
void spfa(int s){
for(int i=1; i<=b[s][0]; i++)
if (dis[b[s][i]] > dis[s]+a[s][b[s][i]]){
dis[b[s][i]] = dis[s]+a[s][b[s][i]];
spfa(b[s][i]);
}
}
int main(){
int x, y, z;
cin >> n >> m;
for(int i=0; i<m; i++){
cin >> x >> y >> z;
if (a[x][y]!=0 && z>a[x][y]) continue;
b[x][0]++; b[x][b[x][0]]=y; a[x][y]=z;
b[y][0]++; b[y][b[y][0]]=x; a[y][x]=z;
}
cin >> s >> t;
for(int i=0; i<=n; i++) dis[i]=99999999;
dis[s]=0;
spfa(s);
if (dis[t]!=99999999) cout << dis[t] << endl;
else cout << -1 << endl;
return 0;
}
spfa優化——前向星優化
星形(star)表示法的思想與鄰接表表示法的思想有一定的相似之處。對每個結點,它也是記錄從該結點出發的所有弧,但它不是采用單向連結清單而是采用一個單一的數組表示。也就是說,在該數組中首先存放從結點1出發的所有弧,然後接着存放從節點2出發的所有孤,依此類推,最後存放從結點n出發的所有孤。對每條弧,要依次存放其起點、終點、權的數值等有關資訊。這實際上相當于對所有弧給出了一個順序和編号,隻是從同一結點出發的弧的順序可以任意排列。此外,為了能夠快速檢索從每個節點出發的所有弧,我們一般還用一個數組記錄每個結點出發的弧的起始位址(即弧的編号)。在這種表示法中,可以快速檢索從每個結點出發的所有弧,這種星形表示法稱為前向星形(forward star)表示法。
例如,在下圖中,仍然假設弧(1,2),(l,3),(2,4),(3,2),(4,3),(4,5),(5,3)和(5,4)上的權分别為8,9,6,4,0,7,6和3。此時該網絡圖可以用前向星形表示法表示如下:
前向星存儲圖:
#include <iostream>
using namespace std;
int first[10005];
struct edge{
int point,next,len;
} e[10005];
void add(int i, int u, int v, int w){
e[i].point = v;
e[i].next = first[u];
e[i].len = w;
first[u] = i;
}
int n,m;
int main(){
int u,v,w;
cin >> n >> m;
for (int i = 1; i <= m; i++){
cin >> u >> v >> w;
add(i,u,v,w);
} //這段是讀入和加入
for (int i = 0; i <= n; i++){
cout << "from " << i << endl;
for (int j = first[i]; j; j = e[j].next) //這就是周遊邊了
cout << "to " << e[j].point << " length= " << e[j].len << endl;
}
}