例 19（最小费用最大流问题）（续例18）由于输油管道的长短不一或地质等原因，
使每条管道上运输费用也不相同，因此，除考虑输油管道的最大流外，还需要考虑输油
管道输送最大流的最小费用。图 8 所示是带有运费的网络，其中第 1 个数字是网络的容
量，第 2 个数字是网络的单位运费。
图8 最小费用最大流问题
解 按照最小费用流的数学规划写出相应的 LINGO 程序如下：
model:
sets:
nodes/s,1,2,3,4,t/:d;
arcs(nodes,nodes)/s 1,s 3,1 2,1 3,2 3,2 t,34,4 2,4 t/:c,u,f;
endsets
data:
d=14 0 0 0 0 -14; !最大流为14;
c=2 8 2 5 1 6 3 4 7;
u=8 7 9 5 2 5 9 6 10;
enddata
min=@sum(arcs:c*f);
@for(nodes(i):@sum(arcs(i,j):f(i,j))-@sum(arcs(j,i):f(j,i))=d(i));
@for(arcs:@bnd(0,f,u));
end
求得最大流的最小费用是 205，而原最大流的费用为 210 单位，原方案并不是最优
的。
类似地，可以利用赋权邻接矩阵编程求得最小费用最大流。LINGO 程序如下：
model:
sets:
nodes/s,1,2,3,4,t/:d;
arcs(nodes,nodes):c,u,f;
endsets
data:
d=14 0 0 0 0 -14;
c=0; u=0;
enddata
calc:
c(1,2)=2;c(1,4)=8;
c(2,3)=2;c(2,4)=5;
c(3,4)=1;c(3,6)=6;
c(4,5)=3;c(5,3)=4;c(5,6)=7;
u(1,2)=8;u(1,4)=7;
u(2,3)=9;u(2,4)=5;
u(3,4)=2;u(3,6)=5;
u(4,5)=9;u(5,3)=6;u(5,6)=10;
endcalc
min=@sum(arcs:c*f);
@for(nodes(i):@sum(nodes(j):f(i,j))-@sum(nodes(j):f(j,i))=d(i));
@for(arcs:@bnd(0,f,u));
end
求最小费用流的一种方法—迭代法
下面我们编写了最小费用最大流函数 mincostmaxflow，其中调用了利用 Floyd 算法
求最短路的函数 floydpath。
求解例 19 具体程序如下（下面的全部程序放在一个文件中）：
function mainexample19
clear;clc;
global M num
c=zeros(6);u=zeros(6);
c(1,2)=2;c(1,4)=8;c(2,3)=2;c(2,4)=5;
c(3,4)=1;c(3,6)=6;c(4,5)=3;c(5,3)=4;c(5,6)=7;
u(1,2)=8;u(1,4)=7;u(2,3)=9;u(2,4)=5;
u(3,4)=2;u(3,6)=5;u(4,5)=9;u(5,3)=6;u(5,6)=10;
num=size(u,1);M=sum(sum(u))*num^2;
[f,val]=mincostmaxflow(u,c)
%求最短路径函数
function path=floydpath(w);
global M num
w=w+((w==0)-eye(num))*M;
p=zeros(num);
for k=1:num
for i=1:num
for j=1:num
if w(i,j)>w(i,k)+w(k,j)
w(i,j)=w(i,k)+w(k,j);
p(i,j)=k;
end
end
end
end
if w(1,num) ==M
path=[];
else
path=zeros(num);
s=1;t=num;m=p(s,t);
while ~isempty(m)
if m(1)
s=[s,m(1)];t=[t,t(1)];t(1)=m(1);
m(1)=[];m=[p(s(1),t(1)),m,p(s(end),t(end))];
else
path(s(1),t(1))=1;s(1)=[];m(1)=[];t(1)=[];
end
end
end
%最小费用最大流函数
function[flow,val]=mincostmaxflow(rongliang,cost,flowvalue);
%第一个参数：容量矩阵；第二个参数：费用矩阵；
%前两个参数必须在不通路处置零
%第三个参数：指定容量值（可以不写，表示求最小费用最大流）
%返回值 flow 为可行流矩阵,val 为最小费用值
global M
flow=zeros(size(rongliang));allflow=sum(flow(1,:));
if nargin<3
flowvalue=M;
end
while allflow<flowvalue
w=(flow<rongliang).*cost-((flow>0).*cost)';
path=floydpath(w);%调用 floydpath 函数
if isempty(path)
val=sum(sum(flow.*cost));
return;
end
theta=min(min(path.*(rongliang-flow)+(path.*(rongliang-flow)==0).*M));
theta=min([min(path'.*flow+(path'.*flow==0).*M),theta]);
flow=flow+(rongliang>0).*(path-path').*theta;
allflow=sum(flow(1,:));
end
val=sum(sum(flow.*cost));