蔬菜运输问题数学建模

1、蔬菜运输问题2012年8月22日摘要本文运用floyd算法求出各蔬菜采购点到每个菜市场的最短运输距离，然后用lingo软件计算蔬菜调运费用及预期短缺损失最小的调运方案，紧接着根据题目要求对算法加以修改得出每个市场短缺率都小于20%的最优调运方案，并求出了最佳的供应改进方案。关键词最短路问题 floyd算法 运输问题一、问题重述光明市是一个人口不到15万人的小城市。根据该市的蔬菜种植情况，分别在花市（A），城乡路口（B）和下塘街（C）设三个收购点，再由各收购点分送到全市的8个菜市场，该市道路情况，各路段距离（单位：100m）及各收购点，菜市场的具体位置见图1，按常年情况，A,B,C三个收购点每天

2、收购量分别为200，170和160（单位：100 kg）,各菜市场的每天需求量及发生供应短缺时带来的损失（元/100kg）见表1.设从收购点至各菜市场蔬菜调运费为1元/(100kg.100m). 7 5 4 8 3 7 A 7 6 B 6 8 5 5 4 7 11 7 4 7 5 6 6 3 5 8 6 6 10 C 10 5 11 图1表1菜市场每天需求（100 kg）短缺损失（元/100kg）7510608805701010010558905808(a) 为该市设计一个从收购点至个菜市场的定点供应方案，使用于蔬菜调运及预期的短缺损失为最小；(b) 若规定各菜市场短缺量一律不超过需求量的20

3、%，重新设计定点供应方案(c) 为满足城市居民的蔬菜供应，光明市的领导规划增加蔬菜种植面积，试问增产的蔬菜每天应分别向A,B,C三个采购点供应多少最经济合理。二、问题分析求总的运费最低，可以先求出各采购点到菜市场的最小运费，由于单位重量运费和距离成正比，题目所给的图1里包含了部分菜市场、中转点以及收购点之间的距离，（a）题可以用求最短路的方法求出各采购点到菜市场的最短路径，乘上单位重量单位距离费用就是单位重量各运输线路的费用，然后用线性方法即可解得相应的最小调运费用及预期短缺损失。第二问规定各菜市场短缺量一律不超过需求量的20%，只需要在上题基础上加上新的限制条件，即可得出新的调运方案。第三问

4、可以在第二问的基础上用灵敏度分析进行求解，也可以建立新的线性问题进行求解。三、模型假设1、各个菜市场、中转点以及收购点都可以作为中转点；2、各个菜市场、中转点以及收购点都可以的最大容纳量为610吨；3、假设只考虑运输费用和短缺费用，不考虑装卸等其它费用；4、假设运输的蔬菜路途中没有损耗；5、忽略从种菜场地到收购点的运输费用。四、符号说明A收购点分送到全市的8个菜市场的供应量分别为a1,b1,c1,d1,e1,f1,g1,h1,B收购点分送到全市的8个菜市场的供应量分别为a2,b2,c2,d2,e2,f2,g2,h2,C收购点分送到全市的8个菜市场的供应量分别为a3,b3,c3,d3,e3,f3

5、,g3,h3,8个菜市场的短缺损失量分别为a,b,c,d,e,f,g,h(单位均为100kg)。五、模型的建立和求解按照问题的分析，首先就要求解各采购点到菜市场的最短距离，在图论里面关于最短路问题比较常用的是Dijkstra算法，Dijkstra算法提供了从网络图中某一点到其他点的最短距离。主要特点是以起始点为中心向外层层扩展，直到扩展到终点为止。但由于它遍历计算的节点很多，所以效率较低，实际问题中往往要求网络中任意两点之间的最短路距离。如果仍然采用Dijkstra算法对各点分别计算，就显得很麻烦。所以就可以使用网络各点之间的矩阵计算法，即Floyd算法。Floyd算法的基本是：从任意节点i到

6、任意节点j的最短路径不外乎2种可能，1是直接从i到j，2是从i经过若干个节点k到j。i到j的最短距离不外乎存在经过i和j之间的k和不经过k两种可能，所以可以令k=1，2，3，.，n(n是城市的数目)，在检查d(i,j)和d(i,k)+d(k,j)的值；在此d(i,k)和d(k,j)分别是目前为止所知道的i到k和k到j的最短距离。因此d(i,k)+d(k,j)就是i到j经过k的最短距离。所以，若有d(i,j)>d(i,k)+d(k,j)，就表示从i出发经过k再到j的距离要比原来的i到j距离短，自然把i到j的d(i,j)重写为d(i,k)+d(k,j)，每当一个k查完了，d(i,j)就是目前

7、的i到j的最短距离。重复这一过程，最后当查完所有的k时，d(i,j)里面存放的就是i到j之间的最短距离了。对于上面的问题，只要能列出它的距离的邻接矩阵，如表2所示。便能用floyd法进行计算了。表2 各点的邻接矩阵图首先对图1中的4个纯中转点进行编号，为（9）（12），并标注在图1中，A、B、C的编号分别为1、14、15，其余点在矩阵中的编号如表1第一行所示。由于数据比较复杂，用普通的计算很困难，所以我们可以用MATLAB软件来编程求解。算法思路：采用标号作业法,每次迭代产生一个永久标号, 从而生长一颗以V0为根的最短路树,在这颗树上每个顶点和根节点之间的路径皆为最短路径。当然此问题也需要表2

8、的各点邻接矩阵。接下来可以得到Dijkstra法的MATLAB语言。M程序如下：function min,path=dijkstra(w,start,terminal)n=size(w,1); label(start)=0; f(start)=start;for i=1:n if i=start label(i)=inf;end, ends(1)=start; u=start;while length(s)<n for i=1:n ins=0; for j=1:length(s) if i=s(j) ins=1; end, end if ins=0 v=i; if label(v)>

9、;(label(u)+w(u,v) label(v)=(label(u)+w(u,v); f(v)=u; end, end, end v1=0; k=inf; for i=1:n ins=0; for j=1:length(s) if i=s(j) ins=1; end, end if ins=0 v=i; if k>label(v) k=label(v); v1=v; end, end, end s(length(s)+1)=v1; u=v1;endmin=label(terminal); path(1)=terminal;i=1; while path(i)=start path(i

10、+1)=f(path(i); i=i+1 ;end path(i)=start;L=length(path);path=path(L:-1:1);图2 Dijkstra算法的MATLAB运行图 如图2，红色框内的是Dijkstra算法的脚本程序，蓝色框内的试运行结果，根据程序显示s点到j点的最短路就是4百米。在程序中显示所走的路线为path=1 2，对应点即为直接从A点到达点。但是由于运用dijkstra编程每次都要修改起始点和终点，所以在大部分情况下效率都不高。由于dijkstra算法较为复杂，所以可用Floyd算法用于求最短路径。Floyd算法思路本质很简单，即用三个for循环的嵌套，代码

11、思路如下：for ( int i = 0; i < 节点个数; +i )for ( int j = 0; j < 节点个数; +j )for ( int k = 0; k < 节点个数; +k )if ( Disik + Diskj < Disij )/ 找到更短路径Disij = Disik + Diskj;但是这里我们要注意循环的嵌套顺序，如果把检查所有节点X放在最内层，那么结果将是不正确的，因为这样便过早的把i到j的最短路径确定下来了，所以正确的应该是这样的： for ( k = 0; k < 节点个数; +k )for ( i = 0; i < 节点

12、个数; +i )for ( j = 0; j < 节点个数; +j )if ( Disik + Diskj < Disij ) 找到更短路径Disij = Disik + Diskj;.那么接下来的问题就是，我们如何找出最短路径.这里需要借助一个辅助数组Path，它是这样使用的：Path(AB)的值如果为P，则表示A节点到B节点的最短路径是A->.->P->B。这样一来，假设我们要找A->B的最短路径，那么就依次查找，假设Path(AB)的值为P，那么接着查找Path(AP)，假设Path(AP)的值为L，那么接着查找Path(AL)，假设Path(AL)的

13、值为A，则查找结束，最短路径为A->L->P->B。那么，如何填充Path的值呢？很简单，当我们发现Dis(AX) + Dis(XB) < Dis(AB)成立时，就要把最短路径改为A->.->X->.->B，而此时，Path(XB)的值是已知的，所以，Path(AB) = Path(XB)。Floyd算法直接在图的带权邻接矩阵中用插入顶点的方法依次递推地构造出n个矩阵D(1), D(2), , D(n), D(n)是图的距离矩阵, 同时引入一个后继点矩阵记录两点间的最短路径。d(i,j) : i到j的距离; path(i,j): i到j的路径上i

14、的后继点;输入带权邻接矩阵a(i,j)。赋初值，对所有i,j, d(i,j)¬a(i,j) , path(i,j)¬j,k=l。然后更新d(i,j) , path(i,j)对所有i,j,若d(i,k)+d(k,j)<d(i,j)则d(i,j)=d(i,k)+d(k,j) , path(i,j)path(i,k) , k =k+1。重复上述步骤直到k=n+1。接下来就是代入具体的数据了，这里的图以及邻接矩阵依旧是修改后的图1及表2。Floyd算法的MATLAB代码如下：M程序function D,path,min1,path1=floyd(a,start,termina

15、l)D=a;n=size(D,1);path=zeros(n,n);for i=1:n for j=1:n if D(i,j)=inf path(i,j)=j;end, end, endfor k=1:n for i=1:n for j=1:n if D(i,k)+D(k,j)<D(i,j) D(i,j)=D(i,k)+D(k,j); path(i,j)=path(i,k);end, end, end,endif nargin=3 min1=D(start,terminal); m(1)=start; i=1; path1= ; while path(m(i),terminal)=ter

16、minal k=i+1; m(k)=path(m(i),terminal); i=i+1; end m(i+1)=terminal; path1=m;end 脚本程序的代码如下：a=0488infinf6infinf7inf4infinfinf407infinfinf5infinfinfinfinfinfinfinf870infinfinfinfinfinf3infinfinf7inf8infinf0inf5infinfinf5inf467infinfinfinfinf0infinfinfinfinfinfinf5infinfinfinfinf5inf0infinfinfinf673inf66

17、5infinfinfinf0infinfinf75infinfinfinfinfinfinfinfinfinf011inf10infinfinf5infinfinfinfinfinfinf110infinfinf6inf107inf35infinfinfinfinf0infinfinf6infinfinfinfinfinf6710infinf0infinfinf84infinf4inf75infinfinfinf0infinfinfinfinfinf653infinf6infinfinf011infinfinf77infinfinfinfinf6infinf110infinfinfinfinf

18、inf6inf510inf8infinfinf0;D, path=floyd(a)运行结果如下：图3 返回矩阵D图4 返回矩阵path D(i,j)表示i到j的最短距离; path(i,j)表示i和j之间的最短路径上顶点i的后继点。根据图3、图4的结果可以很快的知道各点到各个点之间的最短路径，A、（12）、B、C对应1到15这15个网点。例如要找A点到点的最短路，就是从path矩阵寻找。A到，即为1到8，首先找到矩阵中点（1，8）为数字12；再从12出发，找到点（12，8）数字为6；再从6出发，找到点（6，8）为数字15；最后从15出发，找到点（15，8）为数字8。所以最优路线即为从1-12-

19、6-15-8，即从A到（11），（11）到，到C，再从C到，是从A点到点的最短路，其长度可以直接从图3中按（1,8）找得，为22。于是很用以得到图3中红色框内的部分分别就是从A、B、C到各菜市场的最短距离，整理出表3所示的每一百公斤蔬菜从各收购点到各菜市场的运输费用。表3 每一百公斤蔬菜从各收购点到各菜市场的运输费用（元）点购收场市菜A488191162220B14771612162317C20191114615510由于每天三个收购点的总供应量为200+170+160=530（100kg）每天8个菜市场的总需求量为75+60+80+70+100+55+90+80=610（100kg）所以每天

20、的短缺损失量为610-530=80（100kg）（一）对于（a）问题，可以建立以下lindo程序下的线性规划模型：min 4a1+8b1+8c1+19d1+11e1+6f1+22g1+20h1+14a2+7b2+7c2+16d2+12e2+16f2+23g2+17h2+20a3+19b3+11c3+14d3+6e3+15f3+5g3+10h3+10a+8b+5c+10d+10e+8f+5g+8hst 1) a1+b1+c1+d1+e1+f1+g1+h1=2002)a2+b2+c2+d2+e2+f2+g2+h2=1703)a3+b3+c3+d3+e3+f3+g3+h3=1604)a+b+c+d+

21、e+f+g+h=805)a1+a2+a3+a=756)b1+b2+b3+b=607)c1+c2+c3+c=808)d1+d2+d3+d=709)e1+e2+e3+e=10010)f1+f2+f3+f=5511)g1+g2+g3+g=9012)h1+h2+h3+h=80End根据附录1里面的运行结果，我们可以得出各收购点到各菜市场的蔬菜调运量，表4，其最小的蔬菜调运费用及预期短缺损失共为4610元。表4 蔬菜调运量（百公斤）点购收场市菜A750400305500B06040700000C0000700900短缺000000080由于有些收购点到菜市场的最短距离不是直接到达，而是经过其他中转站、菜

22、市场，甚至其他收购点的，因此还要根据图4查出具体的调运线路如下：P（A，1）=A-1，运送量为75百公斤；P（A，3）=A-3，运送量为40百公斤；P（A，5）=A-(11)-5，运送量为30百公斤；P（A，6）=A-6，运送量为55百公斤；P（B，2）=B-2，运送量为60百公斤；P（B，3）=B-3，运送量为40百公斤；P（B，4）=B-（12）-4，运送量为70百公斤；P（C，5）=C-5，运送量为70百公斤；P（C，7）=C-7，运送量为90百公斤；注：P（i,j）为i到j的最短路径因此，按点和点来表示，调运方案还可以这样：表5 点对点的调运方案起点终点（11）（12）A7540553

23、0B604070C7090（11）30（12）70表中单位均为（百公斤），（11）、（12）为中转点，空白处为零，表示无运输需要。（二）显然上一模型得出的结果中菜市场完全没有得到供应，现实生活中是不允许的，那么接下来就解答第二个问题：若规定各菜市场短缺量一律不超过需求量的20%，重新设计定点供应方案。同样的可以用线性模型来求解，于是建立以下lindo程序下的线性规划模型：min 4a1+8b1+8c1+19d1+11e1+6f1+22g1+20h1+14a2+7b2+7c2+16d2+12e2+16f2+23g2+17h2+20a3+19b3+11c3+14d3+6e3+15f3+5g3+10

24、h3+10a+8b+5c+10d+10e+8f+5g+8hst 1) a1+b1+c1+d1+e1+f1+g1+h1=2002)a2+b2+c2+d2+e2+f2+g2+h2=1703)a3+b3+c3+d3+e3+f3+g3+h3=1604)a+b+c+d+e+f+g+h=805)a1+a2+a3+a=756)b1+b2+b3+b=607)c1+c2+c3+c=808)d1+d2+d3+d=709)e1+e2+e3+e=10010)f1+f2+f3+f=5511)g1+g2+g3+g=9012)h1+h2+h3+h=8013)a<1514)b<1215)c<1616)d&l

25、t;1417)e<2018)f<1119)g<1820)h<16End根据附录2里面的运行结果，我们可以得出各收购点到各菜市场的蔬菜调运量，表6，其最小的蔬菜调运费用及预期短缺损失共为4806元。表6 蔬菜调运量（百公斤）点购收场市菜A751000606500B05064560000C00002407264短缺0016141601816同样，根据图4查出具体的调运线路如下：P（A，1）=A-1，运送量为75百公斤；P（A，2）=A-2，运送量为10百公斤；P（A，5）=A-(11)-5，运送量为60百公斤；P（A，6）=A-6，运送量为65百公斤；P（B，2）=B-2，

26、运送量为50百公斤；P（B，3）=B-3，运送量为64百公斤；P（B，4）=B-（12）-4，运送量为56百公斤；P（C，5）=C-5，运送量为24百公斤；P（C，7）=C-7，运送量为72百公斤；P（C，8）=C-8，运送量为64百公斤。因此，按点和点来表示，调运方案如下表7 点对点的调运方案（百公斤）起点终点1112A75100006500600B0506400000056C000024072640011000060000001200056000000（三）要满足城市居民的蔬菜供应同时又要符合经济，就必须是理想的收购量总和等于需求总和，从表7得知，目前的最优运输线路中没有一个收购点是同时作

27、为其他运输线的中转站的，因此，只需把新增加的蔬菜按最优方案加到原来的基础上，而不需要把原来的收购量做减少，可以建立以下lindo程序下的线性规划模型：min 4a1+8b1+8c1+19d1+11e1+6f1+22g1+20h1+14a2+7b2+7c2+16d2+12e2+16f2+23g2+17h2+20a3+19b3+11c3+14d3+6e3+15f3+5g3+10h3st1)a1+a2+a3=752)b1+b2+b3=603)c1+c2+c3=804)d1+d2+d3=705)e1+e2+e3=1006)f1+f2+f3=557)g1+g2+g3=908)h1+h2+h3=809)

28、a1+b1+c1+d1+e1+f1+g1+h1>20010)a2+b2+c2+d2+e2+f2+g2+h2>17011)a3+b3+c3+d3+e3+f3+g3+h3>160end从附录3的运算结果可以得知，增加的80百公斤手工量只需全部分给C收购点，然后重新分配，如表8所示，这时满足了题目要求，同时，蔬菜调运费用及预期短缺损失也是最低的，共4770元。表8 蔬菜调运量（百公斤）点购收场市菜A754000305500B02080700000C00007009080短缺00000000六、模型的检验和改进本文所涉及的最短路问题以及运输供求平衡问题都没有使用常规的图上作业法，而是

29、把算法的思路转化为计算机程序，节省了求解的时间同时也提高了模型的准确度，而且具有较强的可复制性。当然在设计最少运费（最短路）的时候，可以尝试把它和接下来的运输问题用“供求平衡”的方法结合在一起求解，即不用先求出各采购点到菜市场的最短路径，直接按必须供应部分和选择性供应部分用动态规划的方法求解以验证，由于手工的运算量大，在此不作赘述。而第三问的增加供应分配方案，可以从附录2的灵敏度分析中得知改变任何一个的供应量都接改变最终的结果（因为ROW1、2、3的ALLOWABLE INCREASE和ALLOWABLE DECREASE都为零），研究影子价格的意义就不大了。综上所述，本文中所涉及的方法和模型

30、都是合适的。七、模型的推广 本文的解题思路以及所涉及的方法和模型都是准确而且可复制性强的，在解决各种最小费用、最短路线、产销平衡、运输问题时都有较强的参考意义，适当的运用计算机程序解决复杂的计算问题有利于数学使用的推广。八、参考文献【1】运筹学教材编写组，运筹学，清华大学出版社，2011九、附录1、lindo运算结果1（含灵敏度分析） LP OPTIMUM FOUND AT STEP 1 OBJECTIVE FUNCTION VALUE 1) 4610.000 VARIABLE VALUE REDUCED COST A1 75.000000 0.000000 B1 0.000000 0.000

31、000 C1 40.000000 0.000000 D1 0.000000 2.000000 E1 30.000000 0.000000 F1 55.000000 0.000000 G1 0.000000 12.000000 H1 0.000000 5.000000 A2 0.000000 11.000000 B2 60.000000 0.000000 C2 40.000000 0.000000 D2 70.000000 0.000000 E2 0.000000 2.000000 F2 0.000000 11.000000 G2 0.000000 14.000000 H2 0.000000 3

32、.000000 A3 0.000000 21.000000 B3 0.000000 16.000000 C3 0.000000 8.000000 D3 0.000000 2.000000 E3 70.000000 0.000000 F3 0.000000 14.000000 G3 90.000000 0.000000 H3 0.000000 0.000000 A 0.000000 13.000000 B 0.000000 7.000000 C 0.000000 4.000000 D 0.000000 0.000000 E 0.000000 6.000000 F 0.000000 9.00000

33、0 G 0.000000 2.000000 H 80.000000 0.000000 ROW SLACK OR SURPLUS DUAL PRICES 1) 0.000000 -15.000000 2) 0.000000 -14.000000 3) 0.000000 -10.000000 4) 0.000000 -8.000000 5) 0.000000 11.000000 6) 0.000000 7.000000 7) 0.000000 7.000000 8) 0.000000 -2.000000 9) 0.000000 4.000000 10) 0.000000 9.000000 11)

34、0.000000 5.000000 12) 0.000000 0.000000 NO. ITERATIONS= 1 RANGES IN WHICH THE BASIS IS UNCHANGED: OBJ COEFFICIENT RANGES VARIABLE CURRENT ALLOWABLE ALLOWABLE COEF INCREASE DECREASE A1 4.000000 11.000000 INFINITY B1 8.000000 INFINITY 0.000000 C1 8.000000 0.000000 2.000000 D1 19.000000 INFINITY 2.0000

35、00 E1 11.000000 2.000000 0.000000 F1 6.000000 9.000000 INFINITY G1 22.000000 INFINITY 12.000000 H1 20.000000 INFINITY 5.000000 A2 14.000000 INFINITY 11.000000 B2 7.000000 0.000000 INFINITY C2 7.000000 2.000000 0.000000 D2 16.000000 0.000000 2.000000 E2 12.000000 INFINITY 2.000000 F2 16.000000 INFINI

36、TY 11.000000 G2 23.000000 INFINITY 14.000000 H2 17.000000 INFINITY 3.000000 A3 20.000000 INFINITY 21.000000 B3 19.000000 INFINITY 16.000000 C3 11.000000 INFINITY 8.000000 D3 14.000000 INFINITY 2.000000 E3 6.000000 0.000000 2.000000 F3 15.000000 INFINITY 14.000000 G3 5.000000 2.000000 INFINITY H3 10.

37、000000 INFINITY 0.000000 A 10.000000 INFINITY 13.000000 B 8.000000 INFINITY 7.000000 C 5.000000 INFINITY 4.000000 D 10.000000 2.000000 0.000000 E 10.000000 INFINITY 6.000000 F 8.000000 INFINITY 9.000000 G 5.000000 INFINITY 2.000000 H 8.000000 0.000000 INFINITY RIGHTHAND SIDE RANGES ROW CURRENT ALLOW

38、ABLE ALLOWABLE RHS INCREASE DECREASE 1 200.000000 0.000000 0.000000 2 170.000000 0.000000 0.000000 3 160.000000 0.000000 0.000000 4 80.000000 0.000000 0.000000 5 75.000000 0.000000 0.000000 6 60.000000 0.000000 0.000000 7 80.000000 0.000000 0.000000 8 70.000000 0.000000 0.000000 9 100.000000 0.00000

39、0 0.000000 10 55.000000 0.000000 0.000000 11 90.000000 0.000000 0.000000 12 80.000000 0.000000 0.0000002、lindo运算结果2（含灵敏度分析） LP OPTIMUM FOUND AT STEP 20 OBJECTIVE FUNCTION VALUE 1) 4806.000 VARIABLE VALUE REDUCED COST A1 75.000000 0.000000 B1 10.000000 0.000000 C1 0.000000 0.000000 D1 0.000000 2.000000 E1 60.000000 0.000000 F1 55.000000 0.000000 G1 0.000000 12.000000 H1 0.000000 5.000000 A2 0.000000 11.000000 B2 50.000000 0.000000 C2 64.000000 0.000000 D2 56.000000 0.000000 E2 0.000000 2.000000 F2 0.000000 11.000000 G2 0.000000 1