n=5;C=[0 15 16 0 0
0 0 0 13 14
0 11 0 17 0
0 0 0 0 8
0 0 0 0 0]; %弧容量
b=[0 4 1 0 0
0 0 0 6 1
0 2 0 3 0
0 0 0 0 2
0 0 0 0 0]; %弧上單位流量的費(fèi)用
wf=0;wf0=Inf; %wf 表示最大流量, wf0 表示預(yù)定的流量值
for(i=1:n)for(j=1:n)f(i,j)=0;end;end %取初始可行流 f 為零流
while(1)
for(i=1:n)for(j=1:n)if(j~=i)a(i,j)=Inf;end;end;end%構(gòu)造有向賦權(quán)圖
for(i=1:n)for(j=1:n)if(C(i,j)>0&f(i,j)==0)a(i,j)=b(i,j);
elseif(C(i,j)>0&f(i,j)==C(i,j))a(j,i)=-b(i,j);
elseif(C(i,j)>0)a(i,j)=b(i,j);a(j,i)=-b(i,j);end;end;end
for(i=2:n)p(i)=Inf;s(i)=i;end %用 Ford 算法求最短路, 賦初值
for(k=1:n)pd=1; %求有向賦權(quán)圖中 vs 到 vt 的最短路
for(i=2:n)for(j=1:n)if(p(i)>p(j)+a(j,i))p(i)=p(j)+a(j,i);s(i)=j;pd=0;end;end;end
if(pd)break;end;end %求最短路的 Ford 算法結(jié)束
if(p(n)==Inf)break;end %不存在 vs 到 vt 的最短路, 算法終止. 注意在求最小費(fèi)用最大流時(shí)構(gòu)造有
%向賦權(quán)圖中不會(huì)含負(fù)權(quán)回路, 所以不會(huì)出現(xiàn) k=n
dvt=Inf;t=n; %進(jìn)入調(diào)整過(guò)程, dvt 表示調(diào)整量
while(1) %計(jì)算調(diào)整量
if(a(s(t),t)>0)dvtt=C(s(t),t)-f(s(t),t); %前向弧調(diào)整量
elseif(a(s(t),t)<0)dvtt=f(t,s(t));end %后向弧調(diào)整量
if(dvt>dvtt)dvt=dvtt;end
if(s(t)==1)break;end %當(dāng) t 的標(biāo)號(hào)為 vs 時(shí), 終止計(jì)算調(diào)整量
t=s(t);end %繼續(xù)調(diào)整前一段弧上的流 f
pd=0;if(wf+dvt>=wf0)dvt=wf0-wf;pd=1;end%如果最大流量大于或等于預(yù)定的流量值
t=n;while(1) %調(diào)整過(guò)程
if(a(s(t),t)>0)f(s(t),t)=f(s(t),t)+dvt; %前向弧調(diào)整
elseif(a(s(t),t)<0)f(t,s(t))=f(t,s(t))-dvt;end %后向弧調(diào)整
if(s(t)==1)break;end %當(dāng) t 的標(biāo)號(hào)為 vs 時(shí), 終止調(diào)整過(guò)程
t=s(t);end
if(pd)break;end %如果最大流量達(dá)到預(yù)定的流量值
wf=0; for(j=1:n)wf=wf+f(1,j);end;end %計(jì)算最大流量
zwf=0;for(i=1:n)for(j=1:n)zwf=zwf+b(i,j)*f(i,j);end;end %計(jì)算最小費(fèi)用
f %顯示最小費(fèi)用最大流
wf %顯示最小費(fèi)用最大流量
zwf %顯示最小費(fèi)用, 程序結(jié)束
參考資料：圖論算法及其 MATLAB 程序代碼
