基于Busacker-Gowan迭代算法最小费用流代码.zip

function [f wf zwf]=BGf(C,b) % function [f wf zwf]=BGf(C,b) % 计算给定图的最小费用最大流 % C;%弧容量 % b;%弧上单位流量的费用 % f 最小费用最大流矩阵，wf 最大流量 zwf 最小费用 n=size(C,2); wf=0;wf0=Inf; %wf表示最大流量, wf0表示预定的流量值 f=zeros(n,n); %取初始可行流f 为零流 while(1) a=ones(n,n)*inf; for(i=1:n) a(i,i)=0; end%构造有向赋权图 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; %求有向赋权图中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 算法结束 if(p(n)==Inf) break; end %不存在vs 到vt 的最短路, 算法终止。 注意在求最小费用最大流时构造有向赋权图中不会含负权回路, 所以不会出现k=n dvt=Inf;t=n; %进入调整过程, dvt 表示调整量 while(1) %计算调整量 if(a(s(t),t)>0) dvtt=C(s(t),t)-f(s(t),t); %前向弧调整量 elseif(a(s(t),t)<0) dvtt=f(t,s(t)); end %后向弧调整量 if(dvt>dvtt) dvt=dvtt; end if(s(t)==1) break; end %当t 的标号为vs 时, 终止计算调整量 t=s(t); end %继续调整前一段弧上的流f pd=0; if(wf+dvt>=wf0) dvt=wf0-wf;pd=1; end%如果最大流量大于或等于预定的流量值 t=n; while(1) %调整过程 if(a(s(t),t)>0) f(s(t),t)=f(s(t),t)+dvt; %前向弧调整 elseif(a(s(t),t)<0) f(t,s(t))=f(t,s(t))-dvt; end %后向弧调整 if(s(t)==1) break; end %当t 的标号为vs 时, 终止调整过程 t=s(t);end if(pd) break; end%如果最大流量达到预定的流量值 wf=0; for(j=1:n) wf=wf+f(1,j); end end %计算最大流量 zwf=0; for(i=1:n) for(j=1:n) zwf=zwf+b(i,j)*f(i,j); end end%计算最小费用 f; %最小费用最大流