2000年数学建模B题钢管订购和运输

钢管订购和运输

摘要

本文根据问题的条件和要求，建立两个模型，两个模型均为单目标非线性规划模型，并通过

求解这两个模型，完整地解决了问题。

由于铁路运输费用函数具有不可加性，不能直接应用现有的最短路算法来求解铁路和公路交通网中任意两点间最小费用路问题。本文采用了一种分步递推算法，巧妙解决了这一问题。

在单目标非线性规划模型中，将管道铺设分为两个过程。先将钢管从钢管厂运到管道与道路交叉口，再从交叉口铺设到管道线上。这样，总的运输费用就化为两个过程的运输费用之和。本模型是以总费用为目标函数的非线性规划模型，利用Lingo 软件，求出问题一的最优解为1278632万元。

对于问题二通过对模型1的灵敏度分析，确定了S5钢厂的销价的变化对购运计划和总费用的影响最大，确定S1钢厂的生产上限的变化对物运计划和总费用的影响最大。 问题三模型的建立原理和问题一的相同，利用Lingo 软件，求得最优解为1407149万元.

关键词：Floyd算法 单目标非线性规划 灵敏度分析 问题重述

有7个生产厂，可以生产输送天然气主管道的钢管S1,S2,S7。要沿着A1A2A15的主管道铺设, 如题图一所示。图中粗线表示铁路，单细线表示公路，双细线表示要铺设的管道(假设沿管道或者原来有公路，或者建有施工公路)，圆圈表示火车站，每段铁路、公路和管道旁的阿拉伯数字表示里程(单位km)。为方便计，1km主管道钢管称为1单位钢管。 一个钢厂如果承担制造这种钢管，至少需要生产500个单位。钢厂Si在指定期限内能生产该钢管的最大数量为si个单位，钢管出厂销价1单位钢管为pi万元，如下表：siipi1800160280015531000152000160520001556200015073000160 1单位钢管的铁路运价如下表：

里程(km) 运价(万元) 里程(km) ≤300 20 501～600 301～350 23 601～700 351～400 26 701～800 401～450 29 451～500 32 801～900 901～1000 ——仅供参考

运价(万元) 37 44 50

55 60 1000km以上每增加1至100km运价增加5

公路运输费用为1单位钢管每公里0.1万元（不足整公里部分按整公里计算）。 钢管可由铁路、公路运往铺设地点（不只

A1,A2,,A15是运到点，而是管道全线）。

（1）请制定一个主管道钢管的订购和运输计划，使总费用最小（给出总费用)。

（2）请就（1）的模型分析：哪个钢厂钢管的销价的变化对购运计划和总费用影响最大，哪个钢厂钢管的产量的上限的变化对购运计划和总费用的影响最大，并给出相应的数字结果。

（3）如果要铺设的管道不是一条线，而是一个树形图，铁路、公路和管道构成网络，请就这种更一般的情形给出一种解决办法，并对题图二按（1）的要求给出模型和结果。 问题分析 问题一，首先，所有钢管必须运到天然气主管道铺设路线上的节点A1A2A15，然后才能向左或右铺设。必须求出每个钢管厂S1,S2,S7到每个节点A1A2A15的每单位钢管的最小运输费用。 对最小运费的求解，我们采用Floyd算法,先求出铁路网上钢管厂到铁路上任意两点Vi，Vj的最短路线的长度Lij,用matlab求得Lij对应的铁路单位运费Dij；同理用Floyd 算法求出公路网上的任意两点Vj，Vk 的最短公路路线的长度Ljk，结果乘以0.1得到公路运费D1jk。Cikmin(DijD1jk)，j表示所有运输中转点，于是就得到从某钢厂到某铺设点运输单位钢管的最少运输费用。（具体算法及程序见附录） 每个铺设点分别向y,z两个方向展开，通过Lingo编程求出最小铺设费用。运输费用加上购买费用再加上铺设费用就是我们所要求的总费用。 问题二，通过问题一里面Lingo编程运行得出的结果，分析哪个钢厂钢管的销价的变化对购运计划和总费用影响最大，哪个钢厂钢管的产量的上限的变化对购运计划和总费用的影响最大。

问题三，利用同问题一一样的方法，从而可求出某钢厂到某某铺设点运输单位钢管的最少运输费用。（具体算法及程序见附录）

模型的假设与符号说明

1) 基本假设:

○1要铺设的管道侧有公路，可运送所需钢管。

——仅供参考

○2钢管在运输中由铁路运转为公路运时不计中转（换车）费用； ○3所需钢管均由Si(i1,...,7) 钢厂提供； ④假设运送的钢管路途中没有损耗。 2) 符号说明:

Si: 钢厂: 钢厂

SiSi的最大生产能力;

的出厂钢管单位价格(单位: 万元) ;

pid: 公路上一单位钢管的每公里运费(d = 0. 1 万元) ; e: 铁路上一单位钢管的运费(分段函数见表1) ; cijbj: 1 单位钢管从钢厂Si运到Aj的最小费用(单位: 万元) ; : 从Aj 到SiAj1之间的距离(单位: 千米) ; Ajxij: 钢厂运到的钢管数; yj: 运到ZjAjAj地的钢管向左铺设的数目; 地的钢管向右铺设的数目; ：运到1，钢厂Si提供钢管；ti0，钢厂Si不提供钢管； : = W : 所求钢管订购、运输的总费用(单位: 万元) ; 模型的建立与求解 问题一的模型： 针对题图一，我们采用Floyd算法，用matlab编程求出单位钢管从费用，具体数据如下表1：

表1 单位钢管从

SiSi运输到Aj的最小运输运输到

Aj的最小运输费用（单位：万元）

对表1的数据进行分析，我们得到一个非线性规划模型：

目标函数是总费用W , 它包含三项: 钢管出厂总价Q , 运输费P , 及铺设费T. 即 W = Q + P + T

Qpixiji1j1715其中

——仅供参考

，

Pcijxiji1j1715 ，

铺设费T可以如下来确定：

d1d2...dyjd(1yj)yj215Aj开始从左右两个方向铺设，yj与zj单位长钢管的费用为

(1zj)zj2与d

故 Tdj11yjyj1zjzj 22目标函数为： 约束条件为:

500tixijsitij115① 生产能力的: 7) (ti ，(i1,...,0或1) 15) ② 运到Aj的钢管用完: xijyjzj，(j1,...,i1714) ③Aj与Aj-1之间的钢管: zjyj1bj，(j1,...,7,j1,...,15) ④ 变量非负性:xij0,yj0,zj0, (i1,...,⑤ 运到模型一 s.t. 500tixijsitij115Aj的钢管整数: xijN

7) (ti ，(i1,...,0或1) 15) xijyjzj，(j1,...,i1714) zjyj1bj，(j1,..., y1=0 , z15=0 7,j1,...,15) xij0,yj0,zj0, (i1,...,ti=0或1 (i=1,..,7)

d=0.05;

根据模型二编写Lingo程序，程序运行后，得到最优最小费用为W1282142万元。

问题二的模型

——仅供参考

对模型1的灵敏度分析

（1）确定哪个钢厂的销价的变化对购运计划和总费用的影响最大

我们假设该钢厂的销价变化在10%pi万元以内，这是较为合理的，将目标函数的w表示为pi的函数：

因此在销价的变化量相同时，

f越大，则pi的变化对w的变化影响越大。 pif2000单位是最大的，所以S6的销价变化 p6由模型Obj1计算得到的数据可以知道

对购运计算和总费用的影响最大，我们可以通过简单的分析来证明：由于S6提供的数量最 大，销价只要很小变化的，就会引起总费用的很大变化，同时，当价格越来越高，由于P6和x6j互为消长i115的关系，当S6越来越小，它在总需求中占的份额减少，影响减弱，S6下降的速度也将放慢。 除了销价的升高，我们还必须考虑销价的降低，此时应尽量满足提供量最少的点S5，当价格越来越低，由互为消长的关系，S5点的提供量将增加，它在总需求中占的份额增加，影响增强，对于S5上升的速度将放慢。 2）确定哪个钢厂的生产上限的变化对物运计划和总费用的影响最大 由于S1是A1到A8的最优首选，因此若S1与其他Si同时扩大相同的S容量，则S1 会更优，所以推断S1应为影响最大者。由最小费用矩阵C可以知道，Ai(i=1,…,8)所需的 钢管量最好都能由S1提供，则此时S1达到最大需求量，在模型Obj1的条件下，S1为2536 单位，而S1的上限为800单位，考虑到实际钢厂的投入与产出，在很短的时期内生产要达 到原来的3倍，不符合实际意义，所以考虑S1在10%范围内变化。同理对于A0点，最优为S3全部提供，即S3应提供634单位，对于A12,A13,A14,A15,由S6全部提供为最优，即S6应提供1205单位，A11,A10由S5全部提供为最优，即S5应提供796单位。 利用计算机模拟，得到5个供货钢厂分别扩建1%，2%，4%，6%，8%，10%时的成本的增长率，见表 。可以看出，相同的S下A1产生的增长率最大，符合上述分析。一旦工厂扩建范围超过最大需求量，则不再会使目标函数优化，则此时增长率为0。即是上图中S5,S6的情况。而对于S1，一旦S1>2536,则其增长率也为0。（S1的数字结果见表 2 ） △表2 某个Si在变化S的情况下目标函数减小量及减小的比率 △s S1 △z (万元) 872 1744 3488 5232 6976 8720 (%) S2 △z (%) （万元) 328 656 1312 1968 2624 3280 0.025 0.051 0.102 0.153 0.204 0.256 S3 △z (%) （万元) 310 620 1240 1860 2480 3100 0.024 0.048 0.096 0.145 0.193 0.242 S5 △z (万元) 0 0 0 0 0 0 (%) S6 △z (万元) 0 0 0 0 0 0 (%) 1% 2% 4% 6% 8% 10% 0.068 0.136 0.272 0.408 0.4 0.685 0 0 0 0 0 0 0 0 0 0 0 0 问题三的模型

题图二为树形图，采用Floyd算法，用matlab编程求出单位钢管从费用，具体数据如下表2：

——仅供参考

Si运输到

Aj的最小运输

表2 单位钢管从

Si运输到

Aj的最小运输费用 （单位：万元）

Lj由于树形图的出现，则某些管道处会出现多支路。 则模型一中模型的 ，

Rj不再适用，

此时可考虑多增加一些支路变量，并增加约束，在目标函数中增加相应的铺设费。

目标函数： 约束条件：

500tixijsitij121① 生产能力的: ② 运到AjAj (i1,...,7)

(ti0或1)

的钢管用完:xijyjzj (j1,...,21且j9,11,17) i17③ 与Aj114) 之间的钢管:zjyj1bj (j1,...,④ 变量非负性: xij0,yj0,zj0,mj0, (i1,...,7,j1,...,21) ⑤ 运到模型二 s.t. xij0,yj0,zj0,mj0, (i1,...,7,j1,...,21) ti0或1(i=1,..,7) Aj的钢管整数: xijN d=0.05; 根据模型二编写Lingo程序，程序运行后，得到最优最小费用为W1406330万元。

模型优缺点 1. 本文先从简单的角度着手建立模型，采用Floyd算法，简化运输网络。过程严谨，理论性强，逻辑严密，而且易于理解。 2. 在计算最短路径时，我们采用Floyd算法，相比与Dijkstra算法，减少了大量的重复计算，提高了工作效率。

3. 本文大量运用了计算机程序，所有数据均由计算机处理，故误差由计算机精度产生，模型据有良好的稳定性。

参考文献：

[1] 谢金星,薛毅.《优化建模与LINGO/LINGO软件》.北京：清华大学出版社,2005 [2] 宗容,施继红,尉洪,李海燕.《数学实验与数学建模》.云南：云南大学出版社，2009

——仅供参考

[3]陆维新，林皓，陈晓东，《订购与运输钢管的最优方案》.成都：四川大学，6100

附录

附录1

Floyd算法函数在matlab下的M函数文件如下： function [D,path]=floyd(a) n=size(a,1);

D=a;path=zeros(n,n); for i=1:n for j=1:n if D(i,j)~=inf path(i,j)=j; end end end for k=1:n for i=1:n for j=1:n if D(i,k)+D(k,j)——仅供参考

附录2

问题（1）附图1中求最小费用MATLAB程序：

ab=[1 1 2 3 4 5 6 7 8 9 10 11 12 13 15 16 17 18 19 20 20 22 23];

bb=[ 14 15 15 16 19 18 23 24 10 10 11 15 13 14 16 17 19 19 20 21 22 23 24];

w=[20 202 1200 690 690 462 70 30 450 80 1150 1100 306 195 720 520 170 88 160 70 320 160 290]; ab1=[1 2 4 5 6 7 8 9 10 11 14 15 16 17 18 33 34 35]; bb1=[19 20 21 22 23 24 25 26 27 28 29 30 31 31 19 24 31 32]; w1=[3 2 600 10 5 10 12 42 70 10 10 62 30 20 104 31 110 20]; a=sparse(ab,bb,w); a(24,24)=0; a=a+a'; a=full(a); for i=1:24 for j=1:24 if(a(i,j)==0&i~=j) a(i,j)=inf; end end end

[D,path]=floyd(a); a1=sparse(ab1,bb1,w1); a1(35,35)=0; a1=a1+(a1)';

——仅供参考

a1=full(a1); for i=1:35 for j=1:35

if(a1(i,j)==0&i~=j) a1(i,j)=inf; end end end

[D1,path1]=floyd(a1); %距离转换为费用的程序 D1=D1*0.1; %把公路最短距离换算成公路最少费用 for k=1:300 m1(k)={k}; end for k=1:50 m2(k)={300+k}; m3(k)={350+k}; m4(k)={400+k}; m5(k)={450+k}; end for k=1:100

m6(k)={500+k}; m7(k)={600+k};

——仅供参考

m8(k)={700+k}; m9(k)={800+k}; m0(k)={900+k}; end for i=1:24

for j=1:24 %把铁路最短距离换算成铁路最少费用 switch D(i,j) case 0 D(i,j)=0; case m1 D(i,j)=20; case m2 D(i,j)=23; case m3 D(i,j)=26; case m4 D(i,j)=29; case m5 D(i,j)=32; case m6 D(i,j)=37; case m7 D(i,j)=44;

——仅供参考

case m8 D(i,j)=50; case m9 D(i,j)=55; case m0 D(i,j)=60; otherwise D(i,j)=ceil((D(i,j)-1000)/100)*5+60; end end end

%c矩阵表示七个钢管生产厂到十五个铺设节点之间的距离，先把它们都设成20000（任意一个钢管厂到任意一个铺设节点之间的距离不会超过20000），然后用 for 循环求出最小值 c=20000*ones(7,15); for i=1:7 %7个钢管生产厂 for k=18:32 %15个铺设节点 for j=8:24 %7个钢管生产厂和17个中转点，i=1，表示第一个钢管生产厂，j=8，表示第一个中转点 if c(i,k-17)>D(i,j)+D1(j-7,k) c(i,k-17)=D(i,j)+D1(j-7,k);

%对于所有中转点，在铁路网和公路网上的下标相差8 end end

——仅供参考

end end for i=1:7 for k=18:32

if c(i,k-17)>D(i,1)+D1(33,k) c(i,k-17)=D(i,1)+D1(33,k); %33代表第一个钢管生产厂S1点 end if c(i,k-17)>D(i,6)+D1(34,k) c(i,k-17)=D(i,6)+D1(34,k); %34代表第六个钢管生产厂S6点 end if c(i,k-17)>D(i,7)+D1(35,k) c(i,k-17)=D(i,7)+D1(35,k); %35代表第七个钢管生产厂S7点 end end

%因为S1,S6,S7这三个钢管厂有公路直接连接到铺设节点，所以把这三个点单独处理 end

MATLAB程序运行结果： 附录3

问题（1）Lingo程序： model:

——仅供参考

sets:

workplace/1..7/:p,s,t; normdg/1..15/:y,z,b; link(workplace,normdg):c,x; endsets data: d=0.05; s=800 800 1000 2000 2000 2000 3000; b=104,301,750,606,194,205,201,680,480,300,220,210,420,500,0; p=160,155,155,160,155,150,160; c=170.7 160.3 140.2 98.6 38.0 20.5 3.1 21.2 .2 92.0 96.0 106.0 121.2 128.0 142.0 215.7 205.3 190.2 171.6 111.0 95.5 86.0 71.2 114.2 142.0 146.0 156.0 171.2 178.0 192.0 230.7 220.3 200.2 181.6 121.0 105.5 96.0 86.2 48.2 82.0 86.0 96.0 111.2 118.0 132.0 260.7 250.3 235.2 216.6 156.0 140.5 131.0 116.2 84.2 62.0 51.0 61.0 76.2 83.0 97.0 255.7 245.3 225.2 206.6 146.0 130.5 121.0 111.2 79.2 57.0 33.0 51.0 71.2 73.0 87.0 265.7 255.3 235.2 216.6 156.0 140.5 131.0 121.2 84.2 62.0 51.0 45.0 26.2 11.0 28.0 275.7 265.3 245.2 226.6 166.0 150.5 141.0 131.2 99.2 77.0 66.0 56.0 38.2 26.0 2.0; enddata min=w;

w=@sum(link(i,j):(p(i)+c(i,j))*x(i,j))+d*@sum(normdg(j):y(j)^2+y(j)+z(j)^2+z(j)); @for(workplace(i):@sum(normdg(j):x(i,j))>=500*t(i); s(i)*t(i)>=@sum(normdg(j):x(i,j)); @bin(t(i)));

——仅供参考

@for(normdg(j):@sum(workplace(i):x(i,j))=y(j)+z(j)); @for(normdg(j)|j#ne#15:b(j)=y(j)+z(j+1)); z(15)=0;y(1)=0;

@gin(@sum(link(i,j):x(i,j))); end 附录4

问题（3）：附图2中求最小费用MATLAB程序： ab=[1 1 2 3 4 5 6 7 8 9 10 11 12 13 15 16 17 18 19 20 20 22 23]; bb=[ 14 15 15 16 19 18 23 24 10 10 11 15 13 14 16 17 19 19 20 21 22 23 24]; w=[20 202 1200 690 690 462 70 30 450 80 1150 1100 306 195 720 520 170 88 160 70 320 160 290]; ab2=[1 2 4 5 6 7 8 9 10 11 14 15 16 17 18 24 31 32 34 36 36 38]; bb2=[19 20 21 22 23 24 25 26 27 28 29 30 31 31 19 33 34 35 39 37 38 39]; w2=[3 2 600 10 5 10 12 42 70 10 10 62 30 20 104 31 110 20 100 130 190 260]; a=sparse(ab,bb,w); a(24,24)=0; a=a+a'; a=full(a); for i=1:24 for j=1:24

if(a(i,j)==0&i~=j) a(i,j)=inf; end end

——仅供参考

end

[D,path]=floyd(a); a2=sparse(ab2,bb2,w2); a2(39,39)=0; a2=a2+(a2)'; a2=full(a2); for i=1:39 for j=1:39 if(a2(i,j)==0&i~=j) a2(i,j)=inf; end end end for i=17:20 a2(i,i+19)=0; a2(i+19,i)=0; end a2(21,34)=0; a2(34,21)=0; [D2,path2]=floyd(a1);

D2=D2*0.1; %把公路最短距离换算成公路最少费用 for k=1:300 m1(k)={k};

——仅供参考

end for k=1:50

m2(k)={300+k}; m3(k)={350+k}; m4(k)={400+k}; m5(k)={450+k}; end for k=1:100 m6(k)={500+k}; m7(k)={600+k}; m8(k)={700+k}; m9(k)={800+k}; m0(k)={900+k}; end for i=1:24 for j=1:24 %把铁路最短距离换算成铁路最少费用 switch D(i,j) case 0 D(i,j)=0; case m1 D(i,j)=20; case m2 D(i,j)=23;

——仅供参考

case m3 D(i,j)=26; case m4 D(i,j)=29; case m5 D(i,j)=32; case m6 D(i,j)=37; case m7 D(i,j)=44; case m8 D(i,j)=50; case m9 D(i,j)=55; case m0 D(i,j)=60; otherwise D(i,j)=ceil((D(i,j)-1000)/100)*5+60; end end end

%c矩阵表示7个钢管生产厂到21个铺设节点之间的距离，先把它们都设成20000（任意一个钢管厂到任意一个铺设节点之间的距离不会超过20000），然后用 for 循环求出最小值

——仅供参考

c=20000*ones(7,21);

for i=1:7 %7个钢管生产厂 for k=18:39 %21个铺设节点

for j=8:24 %7个钢管生产厂和17个中转点，i=1，表示第一个钢管生产厂，j=8，表示第一个中转点

if (k~=33&k~=34&k~=35&c(i,k-17)>D(i,j)+D2(j-7,k)) c(i,k-17)=D(i,j)+D2(j-7,k); end end end end for i=1:7 for k=18:39 if c(i,k-17)>D(i,1)+D2(33,k) c(i,k-17)=D(i,1)+D2(33,k); %33代表第一个钢管生产厂S1点 end if c(i,k-17)>D(i,6)+D2(34,k) c(i,k-17)=D(i,6)+D2(34,k); %34代表第六个钢管生产厂S6点 end

if c(i,k-17)>D(i,7)+D2(35,k) c(i,k-17)=D(i,7)+D2(35,k);

——仅供参考

%35代表第七个钢管生产厂S7点 end end

%因为S1,S6,S7这三个钢管厂有公路直接连接到铺设节点，所以把这三个点单独处理 end

问题（3）程序运行结果： 附录5

问题（3）Lingo程序： model: sets:

workplace/1..7/:t,p,s; normdg/1..21/:y,z,b,m; link(workplace,normdg):x,c; endsets data:

c=170.7,160.3,140.2,140,38,20.5,3.1,21.2,.2,92,96,106,121.2,128,142, 215.7,205.3,190.2,185,111,95.5,86,71.2,114.2,142,146,156,171.2,178,192, 230.7,220.3,200.2,200,121,105.5,96,86.2,48.2,82,86,96,111.2,118,132, 260.7,250.3,235.2,230,156,140.5,131,116.2,84.2,62,51,61,76.2,83,97, 255.7,245.3,225.2,225,146,130.5,121,111.2,79.2,57,33,51,71.2,73,92, 265.7,255.3,235.2,235,156,140.5,131,121.2,84.2,62,51,45,26.2,11,28, 275.7,265.3,245.2,245,166,150.5,141,131.2,99.2,77,66,56,38.2,22,2; d=0.05;

s=800,800,1000,2000,2000,2000,3000; p=160,155,155,160,155,150,160; b=104,301,750,606,194,205,201,680,300,220,210,420,500,0,0,0,0,0,0,0; enddata min=w;

w=@sum(link(i,j):(p(i)+c(i,j))*x(i,j))+d*@sum(normdg(j)|j#ne#1:(y(j)+y(j)^2)/2)+d*@sum(normdg(j)|j#ne#15:

(z(j)+z(j)^2)/2)+d*(m(9)+m(9)^2)/2+d*(m(11)+m(11)^2)/2+d*(m(17)+m(17)^2)/2; @for(workplace(i):@bnd(500*t(i),@sum(normdg(j):x(i,j)),s(i)*t(i)); @bin(t(i)));

@for(normdg(j)|j#eq#9#or#j#eq#11#or#j#eq#17:@sum(workplace(i):x(i,j))=y(j)+z(j)+m(j));

@for(normdg(j)|j#ne#9#and#j#ne#11#and#j#ne#17:@sum(workplace(i):x(i,j))=y(j)+z(j)); @for(normdg(j)|j#ne#15:z(j)+y(j+1)=b(j)); m(9)+y(16)=42;

——仅供参考

z(19)+y(20)=260; z(20)+z(21)=100; m(11)+m(17)=10; z(17)+y(19)=190; end

z(19)+y(20)=260; z(20)+z(21)=100; m(11)+m(17)=10; y(17)+y(18)=130; z(17)+y(19)=190; end

——仅供参考