道路改造项目中碎石运输的设计

张建方 王 雷 任 祝

( 浙江师范大学 浙江金华 321004)

摘要：本文涉及的是有关道路改造项目中碎石的运输问题，在讨论过程中要涉及到碎石的运输量，以及途中的运输费用。首先建立的是一个优化模型，因此我们便可利用规划来对问题进行求解；再对其目标函数和约束条件进行分析，我们得出的是一组非线性目标函数，所以本题涉及的是一个非线性规划问题。在非线性规划的求解方法上，迄今我们所知道的方法也不多，而且好多方法均具有其约束性，具有普遍性的方法很少。而且一种方法的求解常出现局部解，很容易给我们造成误差。为了避免这种情况，我们需要采用两种方法，比较它们的结果，这样能够给我们提供有力证据。基于这种想法我们采取两种途径解决此问题：1.用LINGO软件解规划，2.采用遗传算法。 关键词：运输费用 优化 非线性规划 遗传算法 文章号：TS200501002

*

一、问题重述（略） 二、模型假设

（1）不计空载时的运输费用。 （2）假设桥的造价接近正无穷，其他价格相对它近似忽略不计。 （3）假设一辆运输车一次铺路方向都仅为单向，不采用左铺和右铺同时进行。 （4）临时道路铺设完马上可以通车，而且运输费与AB间道路的运输费用一致。 （5）河流的宽度足够用于能够在两岸正对面建立两个码头，在计算时不计在河流宽度上的运输费。

三、符号说明 M：道路改造总费用；W：从各采石点运出的碎石量；V1：单位长度的临时道路所占体积；V2：单位长度的公路所占体积；m1：铺设s长度临时道路的运输费用；m2：铺设l长度公路的运输

费用

不考虑码头情况下：

Zi：在AB处的转折点；Z11：过Z1转折点后往左铺路的长度；Z12：过Z1转折点后往右铺路的长度；Z21：过Z2转折点后往左铺路的长度；Z22：过Z2转折点后往右铺路的长度

考虑码头情况下：

Xi：在AB处的转折点；X11：过X1转折点后往左铺路的长度；X12：过X1转折点后往右铺路的长度；X22：过X2转折点后往右铺路的长度；X21：过X2转折点后往左铺路的长度；X31：过X3转折点后往左铺路的长度；X32：过X3转折点后往右铺路的长度

四、问题分析

该问题是一个在一定约束条件下的最优化问题，初步分析题意后可知所求目标函数是个非线性的，难以规划为线性规划问题。由于本题所要涉及的数据变量比较多，所以

*本文获2004年全国大学生电工数学建模竞赛全国二等奖。

很难可以用逐步求精的直接穷举搜索法等简单方法进行求解。此外，题目涉及的费用数据比较大，因此对于最后结果的精度要求不是很高。 由于题目中涉及到坐标轴、抛物线等，再加上数据大，因此在计算方面可能要大量利用到计算机。不过对于非线性规划问题，特别是目标函数复杂的规划，由于理论上也有不完善处，利用计算机软件时不免得到的解具有局部性与不完善性。所以我们将问题转换为非线性规划问题后，首先需要对如何建立模型进行分析，在有码头和无码头的情况下，求出各自的最小值以确定码头是否需建。并且我们要对Lingo得到的解加以检验，避免它只得到局部解。考虑到遗传算法在解决规划问题具有优势，而且可以得到全局解，我们利用C语言编写遗传算法，可以解出另外一个解与Lingo的解进行比较，从而便得到一个相对可靠的解。

五、模型的建立与求解

模型（一）. 根据题意，我们可得： 总费用=运输费+材料成本费（包括码头费） （1）、在无码头的情况下 由于另铺道路进行运输的费用很大，所以在此我们不考虑从某个采石点出发建立两条临时道路，或者途中出现分叉道路。根据题意我们可得：

铺设s长度临时道路的运输费用为：m1=铺设l长度公路的运输费用为： m2参照图（1），我们容易得

∫

=∫0.5*15*20*1000*ldl=75000l

0l0

s

0.1*4*20*1000*sds=4000s2

2

d12=(20-Z1)2+202 d22= (180-Z2)2+572

图（1）

所以，导出目标函数和约束条件如下： min

2

M=7500*(Z112+Z122+Z222+Z212)+4000*(d12+d2)+(Z11+Z12)*1000*7.5*d1*20

+200*1000*7.5*60+(d1+d2)*1000*0.4*60 s.t.

Z11=Z1

Z12+Z21=Z2-Z1

Z1≥50Z2+Z22=200 Z1≤Z2

对上述模型我们采用lingo软件进行求解，其原理很简单，对已知目标函数在约束条件进行局部计算，不到1秒钟时间便可输出结果：

Objective value: 0.22635E+10 Variable Value Reduced Cost Z11 50.00000 0.000000 Z12 70.25767 0.000000 Z21 48.144 -0.1008116E- Z22 31.596 -0.5653175E- Z1 50.00000 0.000000

Z2 168.4031 0.000000

所以从第一行我们便可发现其最小费用为：2.2635*10（元）。（具体程序参见附件[1]）

Z1=50km,Z2=168.4km

9

Z11=50.00000, Z12=70.25767 ,Z21=48.144 ,Z22=31.596 ;

所以Z1位于m4的位置，修建的临时道路上段是沿河切线方向，下段为沿河修建临

时道路，具体方案见图（2）所示。

图（2） （2）、修建码头的情况下

假设需要建立三个码头，因此我们需要在河流上设立C、B、D三个点，作 为码头。

由于两条河流路线近似抛物线，因此可根据题意中给出的m1,m2,m3点得出上游抛物线方程：

y2=−8*(x−50)

根据m5,m6,m7点可导出下游抛物线的解析式为： y2=50*(x−50)

3

2

因此我们可以设B:(50−y2,y2)、C:(50−y1,y1)、D:(50+6y1,y1)，具体见图

2

2

88100

（3）。

图（3）

所以根据总费用公式，我们便可建立目标函数及约束条件：

min

M=4000*(d12+d32+d52+d62)+75000*(X112+X122+X212+X222+X312+X322)+

(X11+X12)*1000*7.5*20*d3+[(X11+X12)*1000*7.5+(d3+d5)*1000*0.4+(X21+X22)*1000*7.5]*(d1*20+d2*6)+[(X21+X22)*1000*7.5+d5*1000*0.4]*d4*6+(X21+X22)*1000*7.5*20*d5+(X31+X32)*1000*7.5*20*d6+200*1000*7.5*60+(d1+d3+d5+d6)*1000*7.5*60

d1=((-30+0.125*y12)2+(20-y1)2)0.5;

d2=0.25*(0.5*y1*(y12+16)0.5+8*log(y1+(y12+16)0.5))-0.25*(0.5*y2*(y22+16)0.5+ 8*log(y2+(y22+16)0.5));d3=((50-0.125*y22-x1)2+y22)0.5;

d4=0.25*(0.5*y2*(y22+16)0.5+8*log(y2+(y22+16)0.5))-0.25*8*log(4)+6/50*(50/6)2*0.5 *log(50/6)-6/50* (0.5*y3*(y32+(50/6)2)0.5+0.5*(50/6)2*log(y3+(y32+(50/6)2)0.5));d5=((50+0.06*y32-x2)2+y32)0.5;d6=((180-x3)2+572)0.5;t.

s.

X11=X1;

X12+X21-X2+X1=0;X22+X31-X3+X2=0;y3≤0;y2≥0;y2≤y1;

X1≤X2;X2≤X3;

X3+X32=200;

通过lingo软件，输出结果为：

Objective value: 0.1762321E+10

Variable Value Reduced Cost Y1 15.48587 0.000000 Y2 15.48587 0.2463340 X1 17.90901 0.000000 Y3 0.000000 0.000000 X2 50.00000 0.000000 X3 171.3401 0.000000 X12 13.49369 0.7579843E- X21 18.59730 0.1723220E- X22 81.97702 -0.3651632E- X31 39.36309 0.1171127E- X32 28.65990 0.000000 X11 17.90901 0.000000

“reduced cost”列出最优单纯形表中判别数所在行的变量的系数，表示当变量有微小变动时，目标函数的变化率，其中基变量的reduced cost 值应为0，对于非基变量Ｘj相应的reduced cost值表示Ｘj增加一个单位（此时假定其他非基变量保持不变）时目标函数减小的量。

具体程序参见附录[2]，通过上表我们不难发现y2=y1,y3=0,所以我们实际只需

建立两个码头便已足够，而且D点位于桥处；最后的最小费用应该还需加上码头的费用，因此最后最小费用为：

1.762321*109+3*105=1.762621*109（元）

通过（1）和（2）进行比较，我们发现造码头后的费用要比不造码头的费用远远要小，而且从程序结果上看，有两个点是重合在一起的，因此我们不需要再考虑设立四个码头的情况。所以为了使总费用最小，需要建立三个码头。码头建立后新的组图参见图（4）。

图（4）

B、C两个码头是建在河的同一点，即两码头建立在河的两正对岸（d2=0）。

从图（4）我们不难得出从采石点S2处运输出去的碎石用于d6、X31、X32的铺设，采石点S1的碎石运往d1、d3、X11、X12、X21、X22的铺设；

WS1=（d1+d3）*V1+（X11+X12+X21+X22）*V2； WS2=d6*V1+（X31+X32）*V2；

代入以上数据得：

WS1=9970（立方米）

WS2=533230（立方米）

各点坐标为：B，C（20.0245，15.48587）、D（50，0）。 模型（二）：

由于在上述模型求解过程中，我们纯采用lingo软件对该模型进行求解，因此对于最后结果的正确性值得商讨。对于非线性规划类的模型，求解方法很多，但均有约束；要得到最正确的答案，其最好的方法就是使用全部搜索法，对于本题由于变量太多，因此我们采用遗传算法。

遗传算法是采用C语言进行编程的一类搜索和优化算法，但它又与传统搜索算法不同,遗传算法从一组随机产生的初始解,称为群体,开始搜索过程。群体中的每个个体是问题的一个解,称为染色体。这些染色体在后续迭代中不断进化,称为遗传。遗传算法主要通过交叉、变异、选择运算实现。交叉或变异运算生成下一代染色体,称为后代。染色体的好坏用适应度来衡量。根据适应度的大小从上一代和后代中选择一定数量的个体,作为下一代群体,再继续进化,这样经过若干代之后,算法收敛于最好的染色体,它很可能就是问题的最优解或次优解。遗传算法中使用适应度这个概念来度量群体中的各个个体的在优化计算中有可能到达最优解的优良程度。度量个体适应度的函数称为适应度函数。适应度函数的定义一般与具体求解问题有关。

其实施过程包括编码、产生初始群体、计算适应度、复制、交换、突变、反复迭代、终止等操作。遗传算法的工作流程图可见下图所示：

基本遗传算法流程图

所以以模型（一）为基础，根据上述目标函数和约束条件，我们可以建立遗传算法的C程序（具体内容参见附录[3]）。 程序的具体设计过程为： （1）编码与约束处理

根据目标函数：min M=f(x1,x2,x3,y1,y2,y3,x12,x21,x22,x31,x32,x11) 建立个体的染色体（即M的可行解）：

X=（x1,x2,x3,y1,y2,y3,x12,x21,x22,x31,x32） 其中各变量受到：

X11=X1;

X12+X21-X2+X1=0; X22+X31-X3+X2=0;

X11=X1;

X12=X2-X21+X2-X1=0;X31=X3-X2-X22

y3≤0;y2≥0;y2≤y1;X1≤X2;X2≤X3;X3+X32=200;

算子修正法处理约束条件：从等式约束消除三个变量得

X32=200-X3 （2）计算适应度与复制

适应度i采用目标函数值i=M

然后采用择优选取法复制出新个体，操作过程：

1）.顺序累计群体中各个体的适应度i，得适应度的累加值

ff

f

Sn：

2）.用

sn=∑fi

Sn除各个体的适应度fi，得相对适应度（m是该群体总数）

Pi=2/m−fi/Sn Pi得

3）.累计

4）.产生[0，1]均匀分布的随机数r;

i，则i个体选入下一代 若i−1

5).反复执行4)得到下一代群体 （3）交换

由于题目的约束条件，我们只运用启发式法进行交换，舍去突变：

gi=∑Pi

g

设

Xi,Xi+1是两父体，Xj为启发式后交换个体，则有：

X=αX+(1−α)X

ii+1

j

α 此处是[0，1]中的随机数

（4）反复（2），（3）直到停止命令 由于遗传算法根据变量多少具有一定的局限性，所以最后的结果与实际数据有点误

差，通过程序运行我们可以得到遗传算法的结果为：2.00936982283*10元（不包括码头费）。具体参数见下图所示：（var(0)—var(10)依次代表

X1=3.5000 X2=61.600 X3=172.400 y1=10.260 y2=3.600

9

y3=-6.400 X1216.199 X2110.401 X22=63.599 X31=47.201 X3227.600

同理我们可以计算出：WS1=952320（立方米）

WS2=584010（立方米）

因此我们可得具体路线图见图（5）所示：

图（5）

各对应点坐标为：C（36.8415，10.260），B（48.38，3.6），D（52.4576，-6.4）

六、模型分析

模型（一）我们采用Lingo软件获得一个最优解，但是对于最终答案的正确与否我们暂且不加以肯定；通过模型（二）我们来进行检验，对于这个模型，我们利用遗传算法，但是遗传算法很受变量个数的，而对于本题问题涉及的变量比较多，而且数据 模型（一）X1=X11=17.91 X2=50 y1=15.48587 y2=15.48587 y3=0

模型（二）X1=X11=35 X2=61.6X3=172.4y1=10.260y2=3.600 y3=−6.400

比较庞大，因此遗传算法最后的结果可能会存在误差。

通过C语言编程后，得到的最后结果比模型（一）略有偏大，所以我们模型（一）最后的结果相对是比较正确的。而遗传算法在这个模型的计算上明显存在误差，具体参数的对照见下表：

模型（一）X12=13.49369 X21=18.59730 X22=81.97702 X31=39.36309 X32=28.65990 模型（二）X12=16.199 X21=10.401 X22=63.599 X31=47.201 X32=27.600

从上表我们可以明显的看出遗传算法最后的结果是存在不合理性的。遗传算法在使用范围上是很广的，但它存在很多局限性，主要是在变量上。由于C语言迭代的次数有限，因此当变量个数多时，它就不能达到全局搜索的功能，不容质疑就有可能使得最后无法搜索到真正的最优解。另外，在约束条件的处理上，当原目标函数的约束条件比较繁时，在遗传算法的编程上就会很难进行处理，也会出现一定误差。

在模型（一）上，经过程序计算，公路AB的长度对最后的最小费用影响最大，因此我们可以得到最小费用随公路长度变化的变化情况：

最小费用/亿元 公路长度/千米

17.62 200

19.00 210

20.47 220

22.04 230

23.72 240

25.51 250

27.41 260

从上表我们也可以看出，最小费用跟公路长度的关系也是处于非线性的，而且是斜率逐渐增大的非线性关系。

在模型（二）中我们给出的是让计算机迭代500次，其实500次对于这样一个多变量的目标函数来说是远远不够的，因此在程序结束时，系统还没有找到真正的最优解，以至于出现最后的误差。但是有点我们必须承认，遗传算法在解决它所能及的问题时，解决的效率是很高的，采用逐步择优法，使你看到它整个得出最优解的过程。对于上述模型也是一样，尽管迭代次数无法满足变量的需求，但是它能清晰的告诉你它计算的始终。应用C程序对模型（二）进行迭代运算（程序略）。从程序结果中可以很清楚的看到，遗传算法是在一个个搜索中不断地接近最优值的算法。而且到最后它能使它所得到的最优解保持稳定。

七、模型总结

本模型是一个涉及到非线性规划的优化模型，在解决非线性规划问题时，最重要的要解决非线性规划的求解方法；我们当然在这方面也下了很大的功夫，对于一个全局的非线性优化模型，而且变量比较多，一般的求解方法是不适用的，当然在软件的使用上

每个软件出来的结果可能会各不相同。在这里我们采用了很普通的一个解规划的软件—Lingo，Lingo软件能解决的变量和约束条件都比较多，而且运行速度比较快。 为了确定最后结果的正确性，我们采用遗传算法对最后结果进行检验，遗传算法是利用C语言进行的一种大范围搜索和优化方法，当然如果能用无穷点搜索法解出来的结果肯定是很准确的。作为遗传算法，它的使用几乎是万能的，但唯一不足之处可能就是在编程上繁琐，对于一个多变量模型，其最后优化的能力也是有限的；所以对于本模型多变量的因素，其结果也肯定是存在误差的。 纵观整个模型，由于知识面的缺乏，我们没能在理论求解方面做得更好，最后的结果尽管比较符合实际，但仍存在不少疑点。在本篇论文上我们的创新之处也就是采用了遗传算法加以对原模型的验证，使得最后的结果更具有说服力，并且我们对两个模型进行了比较更进一步说明了遗传算法的局限性。 对于模型（一）的解决方法，如果在变量少的的情况下，我们可以设置一定步长对其进行一个一个迭代，这是解决非线性规划问题的一个最有效也是最正确的方法；对于模型（二）也是一样，变量多的时候它便会出现一定的误差，本题便是最好的例子；但当然，当数据达到它能处理的范围之内时，遗传算法是一个很具有特色的方法，而且它的应用范围极广。所以在条件允许的情况下，遗传算法是一个很好的数学建模方法，值得提倡。

参考文献： [1]. 赵静等 [2]. 姜启源 [3]. 潘正君等

《数学建模与数学实验》 高等教育出版社 《数学模型》 高等教育出版社 《演化计算》 清华大学出版社