钢管订购和运输的最优规划

1、钢管订购和运输的最优规划摘要本文对钢管订购和运输计划进行了最优化设计。在满足钢厂生产数量约束的条件下，以总费用最小为目标来进行模型建立和求解。本题建模的关键是建立总费用，钢管的订购计划和运输计划三者之间的关系。我们运用图论中的赋权无向图最短路的求解思路，先对图中的结点编号，再以每条路运输费用为权，用Dijkstra算法求图中各钢厂到各管道结点的最小费用路线。将订购的钢管按最小费用路线先运输到管道各结点,然后各结点再向两边铺设。我们运用微元法计算从各管道结点向两边铺设的费用。进而以厂商生产数量等为约束条件，以总费用最小为目标建立了优化模型。为了简化模型求解，根据推出的两个引理，放宽钢厂产量的约束

2、得到新的约束条件，大大减少了模型求解次数。最终用lingo软件对这一模型求解，得到最优订购和运输计划。对于管道是线形时的情况，求得的最小总费用为：1350145万元。在考虑钢厂钢管售价变化和钢管产量上限变化对购运计划和总费用影响方面，我们通过对目标函数求偏导，得到钢厂的销价变化对购运计划和总费用影响最大；通过对影子价格的分析，得到钢厂的产量的上限变化对购运计划和总费用影响最大。对铺设管道是树形图的情形建立了模型，求得最优订购和运输计划，最小费用为1427521万元。本文最后对模型进行了评价，并提出了改进意见。 一 问题重述图一为铺设输送天然气管道的示意图。图中为输送天然气的主管道, 粗线表示铁

3、路，单细线表示公路，双细线表示要铺设的管道(假设沿管道或者原来有公路，或者建有施工公路)，圆圈表示火车站，每段铁路、公路和管道旁的阿拉伯数字表示里程(单位km)。经筛选后可以生产这种主管道钢管的钢厂有。为方便计，1km主管道钢管称为1单位钢管。一个钢厂如果承担制造这种钢管，至少需要生产500个单位。钢厂在指定期限内能生产该钢管的最大数量为个单位，钢管出厂销价1单位钢管为万元。 运价分为铁路运价和公路运价两种。对于铁路运价，1单位的钢管在不同的铁路路程范围内对应不同的运价，1000km以上每增加1至100km运价增加5万元；公路运输费用为1单位钢管每公里0.1万元（不足整公里部分按整公里计算）。

4、钢管可由铁路、公路运往铺设地点（不只是运到点，而是管道全线）。需要制定一个主管道钢管的订购和运输计划，使总费用最小（给出总费用）；并就(1)的模型分析哪个钢厂钢管的销价的变化对购运计划和总费用影响最大，哪个钢厂钢管的产量的上限的变化对购运计划和总费用的影响最大，并给出相应的数字结果；同时需进一步讨论要铺设的管道是一个树形图的情况，就这种更一般的情形给出一种解决办法，并对图二按（1）的要求给出模型和结果。1） 铁路运费求解程序：traincost.m2） 最小费用路线的费用求解程序：Dij.mDwij2.mDwij2final.m3） 最小费用路线的路线求解程序：deallamda.m4） 显示

5、结果程序： ShowD.mshowlj.m5） 数据处理程序：dealdata.m6） lingo求解程序traincost.mfunction Tc=traincost(x)L=300 350 400 450 500 600 700 800 900 1000;C=20 23 26 29 32 37 44 50 55 60;if(x>L(end) Tc=C(end)+5*ceil(x-1000)/100);else L(end+1)=x; L=sort(L); xsub=find(L=x); if(xsub(1)=1) Tc=C(1); else if(L(xsub(1)-1)=L(xs

6、ub) Tc=C(xsub(1)-1); else Tc=C(xsub(1); end endendDij.mfunction C,LAMDA=Dij(w)C=ones(7,14)*Inf;w1=w+w'LAMDA=zeros(7,38);for D=1:7 T=ones(1,38)*Inf; P=zeros(1,38); lamda=zeros(1,38); S=zeros(1,38); k=D; S(D)=k; i=1; lamda(D)=-1; while(1) if(length(find(S=0)=0)%判断S是否都已经解决 break; end Lsub=0; sub=0;

7、 sub=find(w1(k,:)=0);%查找k行所有的相临的节点 if(length(find(S(sub)=0)=0)%对相临节点还没便历过的 for j=1:length(sub) %对每一个相临节点遍历， if(S(sub(j)=0)%只有没检查时才检查，如果已经检查过了就不检查 if(sub(j)<25&k<25)%开始计算费用，判断当前节点与前一节点为铁路时，进行铁路运费的计算 Lsub(1)=sub(j); Lsub(2)=k; v=3; ktemp=k; while(lamda(ktemp)<25&lamda(ktemp)=-1)%统计到不是

8、铁路的节点为止 Lsub(v)=lamda(ktemp); ktemp=Lsub(v); v=v+1; end costsum=0;%记录路线中的铁路费用 for ii=1:length(Lsub)-1%所有的铁路段相加，求出铁路长度的总和 costsum=w1(Lsub(ii),Lsub(ii+1)+costsum; end costsum=traincost(costsum);%调用traincost求出铁路费 totalcost=costsum+P(Lsub(end);%记录从起点开始到这点的总费用 else totalcost=P(k)+w1(k,sub(j);%如果是公路，只需要相加

9、即可 end%end if sub if(T(sub(j)>totalcost)%如果T标号中的对应点比求出的总费用大，则替换T标号中的内容 T(sub(j)=totalcost; lamda(sub(j)=k; end else subtemp(j)= 0; end%end if S end%end forj elseif(length(find(S(sub)=0)=0)%没有全部找到,没有完成解题的要求，但是相临的点都已经找到，那么用前一次的起点 %k=karray(i-1); Sv=S(find(S=0); %得到S中所有不等于0，主要为了查找另一个起点 flag=0;%找点的标志

10、，0未找到点,1找到 for Svi=1:length(Sv)%遍历每一个起点，查找还有可以继续进行计算的相临的点作为起点 sub2=find(w1(Sv(Svi),:)=0);%查找已经标号的点的相临的点。 for sub2i=1:length(sub2)%看是否他的每一个相临点都已经找过了 ftemp=find(w1(sub2(sub2i),:)=0);%相临点与哪些点相连% isnotdian=;l& if(length(find(S=sub2(sub2i)=0&sub2(sub2i)%这里排除了一个顶点不在有延伸的点，且相连的一点以遍历 if(length(ftemp)

11、=1&ftemp=Sv(Svi) k=sub2(sub2i); %最终在这里找到k点 flag=1; break; end end end if(flag=1) %如果找到则跳出这个循环 break; end end end%end if length % Tv=min(T(subtemp(find(subtemp=0); Tv,Tvsub=min(T);% if(Tv<Inf&length(find(subtemp=0)=0) if(Tv<Inf) T(Tvsub)=Inf; P(Tvsub)=Tv;%相应的最小值存入相应的P标号当中 S(Tvsub)=Tvsub

12、;%同时把这个找到的最小值的相临的点放入S中 k=Tvsub;%把找到的相临的点作为下一个起始点 elseif(length(find(S(sub)=0)=0)%没有全部找到,没有完成解题的要求，但是相临的点都已经找到，那么用前一次的起点 %k=karray(i-1); Sv=S(find(S=0); %得到S中所有不等于0，主要为了查找另一个起点 flag=0;%找点的标志，0未找到点,1找到 for Svi=1:length(Sv)%遍历每一个起点，查找还有可以继续进行计算的相临的点作为起点 sub2=find(w1(Sv(Svi),:)=0);%查找已经标号的点的相临的点。 for su

13、b2i=1:length(sub2)%看是否他的每一个相临点都已经找过了 ftemp=find(w1(sub2(sub2i),:)=0);%相临点与哪些点相连% isnotdian=;l& if(length(find(S=sub2(sub2i)=0&sub2(sub2i)%这里排除了一个顶点不在有延伸的点，且相连的一点以遍历 if(length(ftemp)=1&ftemp=Sv(Svi) continue; end k=sub2(sub2i); %最终在这里找到k点 flag=1; break; end end if(flag=1) %如果找到则跳出这个循环 bre

14、ak; end end else break;% out while所有点都找到的时候就结束此问题的一个工厂起点的搜索 end%end if Tv i=i+1;%i递增继续 end%end while C(D,:)=P(25:38);% lamdalength=size(lamda,2);% LAMDA(D,lamdalength)=0;if(size(lamda,2)=38) lamda(1,38)=0;end LAMDA(D,:)=lamda;end%end for DDwij2.mfunction C,LAMDA,wtemp=Dwij2(w)%计算从S1,S2,S3,S4,.S7到各个铁

15、路节点的最小费用得到 7*24的矩阵C=ones(7,24)*Inf;wtemp=w;wtemp(1:24,1:24)=0;w=w(1:24,1:24);%只需要所有的铁路的接点的信息w1=w+w'LAMDA=zeros(7,24);for D=1:7 T=ones(1,24)*Inf; P=zeros(1,24); lamda=zeros(1,24); S=zeros(1,24); k=D; S(D)=k; i=1; lamda(D)=-1; while(1) if(length(find(S=0)=0)%判断S是否都已经解决 break; end sub=0; Lsub=0; su

16、b=find(w1(k,:)=0);%查找k行所有的相临的节点 if(length(find(S(sub)=0)=0)%对相临节点还没便历过的 for j=1:length(sub) %对每一个相临节点遍历， if(S(sub(j)=0)%只有没检查时才检查，如果已经检查过了就不检查 %开始计算费用，判断当前节点与前一节点为铁路时，进行铁路运费的计算 Lsub(1)=sub(j); Lsub(2)=k; v=3; ktemp=k; while(lamda(ktemp)=-1)%统计到起点为止 Lsub(v)=lamda(ktemp); ktemp=Lsub(v); v=v+1; end cos

17、tsum=0;%记录路线中的铁路费用 totalcost=0; for ii=1:length(Lsub)-1%所有的相加，求出 总费用 costsum=w1(Lsub(ii),Lsub(ii+1)+costsum; end totalcost=traincost(costsum);%调用traincost求出铁路费%记录从起点开始到这点的总费用 if(T(sub(j)>totalcost)%如果T标号中的对应点比求出的总费用大，则替换T标号中的内容 T(sub(j)=totalcost; lamda(sub(j)=k; end end%end if S end%end forj end

18、%end if length Tv,Tvsub=min(T); if(Tv<Inf) T(Tvsub)=Inf; P(Tvsub)=Tv;%相应的最小值存入相应的P标号当中 S(Tvsub)=Tvsub;%同时把这个找到的最小值的相临的点放入S中 k=Tvsub;%把找到的相临的点作为下一个起始点 else break;% out while所有点都找到的时候就结束此问题的一个工厂起点的搜索 end%end if Tv i=i+1;%i递增继续 end%end while C(D,:)=P; if(size(lamda,2)=24) lamda(1,24)=0; end LAMDA(D,

19、:)=lamda; end%end for D %的得出总的矩阵 wtemp(1:7,1:24)=C; wtemp=wtemp+wtemp' wtemp(1:7,1:7)=C(1:7,1:7);Dwij2finalfunction C,LAMDA=Dwij2final(w1)%here w1 is waiting for Dwij2's return value wtemp as nargin w1C=ones(7,14)*Inf;LAMDA=zeros(7,38);for D=1:7 T=ones(1,38)*Inf; P=zeros(1,38); lamda=zeros(1

20、,38); S=zeros(1,38); k=D; S(D)=k; i=1; lamda(D)=-1; while(1) if(length(find(S=0)=0)%判断S是否都已经解决 break; end sub=0; Lsub=0; sub=find(w1(k,:)=0);%查找k行所有的相临的节点，得到点的标号； % k,sub,S(sub),S(sub)=0,find(S(sub)=0) if(isempty(sub) if(length(find(S(sub)=0)=0)%对相临节点还没便历过的 for j=1:length(sub) %对每一个相临节点遍历， if(S(sub(

21、j)=0)%只有没检查时才检查，如果已经检查过了就不检查 %开始计算费用，进行运费的计算,直接相加就可以了 Lsub(1)=sub(j); Lsub(2)=k; v=3; ktemp=k; while(lamda(ktemp)=-1)%统计到起点为止 Lsub(v)=lamda(ktemp); ktemp=Lsub(v); v=v+1; end costsum=0;%记录路线中的铁路费用 totalcost=0; for ii=1:length(Lsub)-1%所有的相加，求出 总费用 costsum=w1(Lsub(ii),Lsub(ii+1)+costsum; end totalcost=

22、costsum;% if(T(sub(j)>totalcost)%如果T标号中的对应点比求出的总费用大，则替换T标号中的内容 T(sub(j)=totalcost; lamda(sub(j)=k; end end%end if S end%end forj end%end if length end %end if sub=null Tv,Tvsub=min(T); if(Tv<Inf) T(Tvsub)=Inf; P(Tvsub)=Tv;%相应的最小值存入相应的P标号当中 S(Tvsub)=Tvsub;%同时把这个找到的最小值的相临的点放入S中 k=Tvsub;%把找到的相临的点

23、作为下一个起始点 else break;% out while所有点都找到的时候就结束此问题的一个工厂起点的搜索 end%end if Tv i=i+1;%i递增继续 end%end while C(D,:)=P(25:38);% lamdalength=size(lamda,2);% LAMDA(D,lamdalength)=0;if(size(lamda,2)=38) lamda(1,38)=0;end LAMDA(D,:)=lamda;end%end for Ddeallamda.mfunction Passby,time2,Ptemp = deallamda( LAMDA )%UNTI

24、TLED1 Summary of this function goes here% Detailed explanation goes here%在一列中分别对应 A1 ,A2,.,A15。其中A1未编号。程序要根据编号来进行考虑times=10;%存放第二维的大小time2=0;Passby=zeros(14,times,7);%considerate 10 times is enoughfor i=1:7 P=zeros(1,times); lamda=LAMDA(i,:); temp=0;for j=25:38 %compute every point of S ,that is,How

25、 they get to the A points; temp=lamda(j); k=1; while(temp=-1)%进行搜索，直到搜到起点起点的lamda为1 if(k>times) time2=k; end P(k)=temp; temp=lamda(temp); k=k+1; end if(time2>times) %动态的设置维数 Passby(14,time2,7)=0; times=time2; end Passby(j-24,:,i)=P; P=zeros(1,times);endendPtemp=zeros(15,times+1,7);tempw=ones(1

26、,1,7)*25;Ptemp(1,:,:)=tempw,Passby(1,:,:);Passby(14,times+1,7)=0;Ptemp(2:15,:,:)=Passby(:,:,:);ShowD.mfunction ShowD( Ptemp )%UNTITLED1 Summary of this function goes here% Detailed explanation goes hereu,v,w=size(Ptemp);for wi=1:w P=Ptemp(:,:,wi); disp(sprintf('S%d开始:n',wi); for ui=1:u P2=P(

27、ui,:); subtemp=find(P2=wi); P2(subtemp(1)+1:end)=0; P2(find(P2=0)=; P2=fliplr(P2); P3=num2str(P2); P4=strrep(P3,' ',' >'); disp(sprintf('t到A%d :',ui); disp(sprintf('t%s',P4); endendshowlj.mload wx,y=Dwij(w);u,v,ww=deallamda(y);ShowD(ww);Dealdata.mfunction C = deald

28、ata( D )%compute the cost of the unit fee ,that the A will pay for the thing they%get% Detailed explanation goes hereP=160;155;155;160;155;150;160;C=ones(7,15);C(:,1)=D(:,1)+10.4;C(:,2:15)=D;C=C+P(:,ones(1,15);m,n=size(C);fid=fopen('d:C2.doc','w+');for i=1:mfor j=1:nif (j=n)fprintf(f

29、id,'%4.1ft',C(i,j);elsefprintf(fid,'%4.1fn',C(i,j);endendendfclose(fid);程序使用说明：两种方法进行 1，运行deal3,输入D,得到单价的矩阵7*15;运行showlj,得到运输的路线； 2, 运行deal2,输入D,得到单价的矩阵7*15;deal2.mload w;C,LAMDA,wtemp=Dwij2(w);C,LAMDA=Dwij2final(wtemp); D = dealdata( C );deal3.mload w;C,LAMDA=Dij(w); D = dealdata( C );lingo 求解程序 model:sets: numberi/1.7/:s,w; numberj/1.15/:m1,m2,m,l; numberc(numberi,numberj):c,t;endsetsdata: l=0,104,301,750,606,194,205,201,680,480,300,220,210,420,500; c=330.7320.3300.2258.6198.0180.5163.1181.2224.2252.0256.0266.0281.2288.0302.0370.7360.3345.2326.