第二讲数学建模竞赛中的部分优化题目

数学建模竞赛中的部分问题描述0（题数学建模竞赛中的部分问题描述0（题(i1,2,,6)0优化问题

2.1一个飞行管理问题

2.1.1

1995年全国大学生数学建模竞赛中的A题（“一个飞行管理问题”）。在约10000米高空的某边长为160km的正方形区域内，经常有若干架飞机做水平飞行，区域内每架飞机的位置和速度向量均由计算机记录其数据，以便进行飞行管理，当一架欲进入该区域的飞机到达区域边缘时，记录其数据后，要立即计算并判断是否会与区域内的飞机发生碰撞。现假定条件如下：（1）不碰撞的标准为任意两架飞机的距离大于8km；（2）飞机飞行方向角调整幅度不应超过30°；（3）所有飞机飞行速度均为800km／h；（4）进入该区域的飞机在到达该区域边缘时，与区域内飞机的距离应在60km以上；（5）最多需考虑6架飞机；（6）不必考虑飞机离开此区域后的情况。请你对这个避免碰撞的飞行管理问题家里数学模型，列出计算步骤，对以下数据进行计算（方向角误差不超过0.01°），要求飞机飞行方向角调整的幅度尽量小。该区域四个定点的坐标为（0，0）、（160，0）、（160，160）、（0，160）。记录数据见表2－1。表2－1飞机位置和方向角记录数据飞机横坐纵坐方向角飞机横坐纵坐方向角／编号标标（°）编号标x标y（°）11501402434145501592858523651301502303155155220．5新进0052入说明：方向角指飞行方向与x轴正向的夹角。试根据实际应用背景对你的模型进行评价和推广。

2．1．2模型一及求解

模型建立这个问题显然是一个优化问题。设第i架飞机在调整时的方向角为i目中已给出），调整后的方向角为，题目中就是要iii求飞机飞行方向角调整的幅度尽量小，因此有化的目的函数可以是：

||2..v（＝800km／h），以当前时刻为

(x0,y0)，则Ttir(t)，这时在区域内飞机不相撞的约束条件就变成了（5）||2..v（＝800km／h），以当前时刻为

(x0,y0)，则Ttir(t)，这时在区域内飞机不相撞的约束条件就变成了（5）（6）

ijij2ijij22（1）（2）（3）（7）（8）（9）ii1为了建立这个问题的优化模型，只需要明确约束条件就可以了。一个简单的约束是飞机飞行方向角调整的幅度不应超过30°，即

||30题目中要求进入该区域的飞机在到达该区域边缘时，与区域内的飞机的距离应在60km以上。这个条件是个初始条件，很容易验证目前所给的数据是满足的，因此本模型中可以不予考虑。剩下的关键是要慢足题目中描述的任意两架位于该区域内的飞机的距离应该大于8km。但这个问题的难点在于飞机是动态的，这个约束不好直接描述，为此我们首先需要描述每架飞机的飞行轨迹。

记飞机飞行速率为0时刻。设第i架飞机在调整时的位置坐标为（已知条件），t时刻的位置坐标为ii

(xt,yt)ii

xtx0vtcos,yty0vtsin.iiiiii如果要严格表示两架位于该区域内的飞机的距离应大于8km，则需要考虑每架飞机在区域内的飞行时间的长度。记为第架飞机飞出区域i的时间，即

Targmin{t0:x0vtcos0或160，（4）iii或y0vtsin0或160｝ii

记时刻第架飞机与第架飞机的距离为，并记ij

f(t)[r(t)]264ijij

f(t)[r(t)]2640(0tT)ijijij其中

Tmin{T,T}.ijij此外，经过计算可以得到

f(t)(x0vtcosx0vtcos)2ijiijj(y0vtsiny0vtsin)264z2bzciijjijijij

z2vtsin,ij

b2[(x0x0)sin(y0y0)cos],ijijij

（10）f(t)zijijijij。注意到）则此时约束条件（5）一定成立，所以0ijij（11）0ijijijijijijij（12）f(t)6（14）（10）f(t)zijijijij。注意到）则此时约束条件（5）一定成立，所以0ijij（11）0ijijijijijijij（12）f(t)6（14）（15）（16）b/2,tb/4vsintf(t)且t，只要在右端点的函数值非负即可，即t*Tf(t)f(t)b/4c0||2;即（记为）时，函数取最小值2*，只需要的最小值（13）ij*T*ijijij

所以是一个关于的二次函数，表示的是一条开口向上的抛物线。ij

当ij

b/4cf(0)c0（初始时刻不相撞），如果b/20（即ijijijijij

b0ij如果bij

f(T)0;ijij如果b且0ij即可，即

b24c0.ijij

实际上，约束（11）表示的是在右端点的函数值非负，这个约ij束在（12）的条件下也是自然成立的，所以可以是对约束（11）不再附加且的条件。于是我们的模型就是

minii1

s.t.||230,i

f(T)0;ijij

b24c0.(b0,0t*T)ijijijijij

模型求解

上面这是一个非线性规划模型，虽然是严格满足题目要求的模型，但得到的模型逻辑关系比较复杂，约束（16）是在一定条件下才成立的约束，而且其中的计算式（4）也含有相当复杂的关系式，使用LINGO软件不太容易将模型很方便的输入，因为逻辑处理不是LINGO的优势所在。即使想办法把这个模型输入到LINGO，也不一定能求出好的解（笔者尝试过，但是LINGO运行时有时会出现系统内部错误，可能是系统有问题，无法继续求解）。而且，在实时飞行调度中显然需要快速求解，所以下面

1602/8000.21602/8000.220.283(h)T2，任何一架飞机在所考虑的区域内这个模型麻烦之处就在于，要求严格表示两架飞机的飞行距离应大于8km，所以需要考虑每架飞机在区域内的飞行时间的长度，比较繁琐。

注意到区域对角线的长度只有

停留的时间不会超过Tmax160。因此这里我们简化一下问题；不再单独考虑每架飞机在区域内停留的时间，而是以最大时间T（这是已经是一个常数）代替之，此时所有T，这实际maxijmax上强化了问题的要求，即考虑了有些飞机可能已经飞出区域，但仍不允许两架飞机的距离小于8km。

程序：

MODEL:TITLE飞行管理问题的非线性规划模型;SETS:Plane/1..6/:x0,y0,cita0,cita1,d_cita;!cita0表示初始角度，cita1为调整后的角度,d_cita为调整的角度;link(plane,plane)|&1#LT#&2:b,c;ENDSETSDATA:x0y0cita0=1501402438585236150155220.5145501591301502300052;max_cita=30;T_max=0.283;V=800;ENDDATAINIT:d_cita=000000;ENDINIT@for(plane:cita1-cita0=d_cita);@for(link(i,j):b(i,j)=-2*(x0(i)-x0(j))*@sin((cita1(i)+cita1(j))*3.14159265/360)+2*(y0(i)-y0(j))*@cos((cita1(i)+cita1(j))*3.14159265/360);c(i,j)=(x0(i)-x0(j))^2+(y0(i)-y0(j))^2-64;);!避免碰撞的条件;

i66|

i66||||2@for(link(i,j):[Right](2*V*T_max*@sin((cita1(i)-cita1(j))*3.14159265/360))^2+b(i,j)*(2*V*T_max*@sin((cita1(i)-cita1(j))*3.14159265/360))+c(i,j)>0);!最小点非负;@for(link(i,j):[Minimum]@if(b(i,j)#lt#0#and#

-b(i,j)/4/V/@sin((cita1(i)-cita1(j))*3.14159265/360)#gt#0#and#

-b(i,j)/4/V/@sin((cita1(i)-cita1(j))*3.14159265/360)#lt#T_max,b(i,j)^2-4*c(i,j),-1)<0);!@for(link(i,j):@if(b(i,j)#lt#0,b(i,j)^2-4*c(i,j),-1)<0);@for(link:@free(b));!调整角度上下限，单位为角度;@for(plane:@bnd(-max_cita,d_cita,max_cita));[obj]MIN=@SUM(plane:(d_cita)^2);![obj]MIN=@SUM(plane:@abs(d_cita));END

运行这个程序，结果得到的是一个局部极小点，调整角度较大。能找到更好的解吗？如果不用全局求解程序，通常很难得到稍大规模的非线性规划问题的全局最优解。所以我们启动LINGO全局求解程序求解这个模型（过程省略），可以得到全局最优解。可以看出，在0。01度的误差要求下，需要调整第3、4、6三架飞机的角度，分别调整2.06度，－0.5度,1.57度，调整量的平方和为6.95。其实，使用全局变量求解程序，通常也不一定要等到全局最优解，而是观察求解状态窗口，看到一个较好的当前解（或当前最好解在较长时间内不发生变化）时，就可以中止程序，用当前最好的局部最优解作为最后结果。例如对于本例，LINGO求出最优解大约需要一分，而实际上5秒内LINGO得到了与全局最优解类似的解。此外，上面的模型还可以进一步简化，例如可以假设要求飞机永远不相撞，即认为为无穷大，这时显然约束条件（15）是多余的，而且约束（16）中只需要的条件就可以了。也就是说上面程序中的对应部分（约束和可以改写为更简单的形式：！右端点非负，不再需要；！最小点非负，需化为以下形式：实际计算结果显示，此时得到的结果与前面计算的结果几乎没有差别。

备注：优化的目标函数除了外，也可以假定为或iii1i1

max||等，用LINGO求解的过程的是完全类似的，计算结果略有差异，1i6

AAA2S,S,S2Ssp1800160≤30020501～600371572800155301～35023601～7004431000155351～40026701～8005042000160401～45029801～900AAA2S,S,S2Ssp1800160≤30020501～600371572800155301～35023601～7004431000155351～40026701～8005042000160401～45029801～9005552000155451～50032901～1000606200015073000160机的数量尽量少，这种想法在实际中也不能说没有道理，但与题目的要求不符，而且解题难度并没有减小，意义似乎不大。在实际调度中，由于计算上面的调度方案，需要时间将调度信息告知飞机驾驶员并做出调整方向角的操作也需要时间，因此如果考虑一定的反应滞后时间应该是比较合理的。也就是说，如果反应时间是10秒，则计算中应采用飞机沿当前方向角飞行10秒以后的位置作为计算的基础。2.2钢管订购和运输

2.2.1问题描述

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

i

si

pi

1单位钢管的铁路运价如下表：

里程(km)运价(万元)

里程(km)运价(万元)

1000km以上每增加1至100km运价增加5万元。公路运输费用为1单位钢管每公里0.1万元（不足整公里部分按整公里计算）。

A,A,230S46903070170520462S14803151080750A2,A15S7160S611088S57042A910194A6606A4301，而32070A156210A1310A10201205A5A316020500420220300A8A7A,A,230S46903070170520462S14803151080750A2,A15S7160S611088S57042A910194A6606A4301，而32070A156210A1310A10201205A5A316020500420220300A8A7图一20A14210A12A116801是管道全线）。（1）请制定一个主管道钢管的订购和运输计划，使总费用最小（给出总费用)。（2）请就（1）的模型分析：哪个钢厂钢管的销价的变化对购运计划和总费用影响最大，哪个钢厂钢管的产量的上限的变化对购运计划和总费用的影响最大，并给出相应的数字结果。（3）如果要铺设的管道不是一条线，而是一个树形图，铁路、公路和管道构成网络，请就这种更一般的情形给出一种解决办法，并对图二按（1）的要求给出模型和结果。

290

S3S2

6901200720

202

110020121953061150

600450

23104

A1

17S1480A113168051080

750

301A3A2A,A,,A2S,S,7（记为j=1，2，„„，15）的每单位钢管的最小运费Cij（不21601001308846242A910194A6

415S1532070190S517S1480A113168051080

750

301A3A2A,A,,A2S,S,7（记为j=1，2，„„，15）的每单位钢管的最小运费Cij（不21601001308846242A910194A6

415S1532070190S57010300201205

606（记为i=1，2，„„，7）到每个结点160A20260(AA191010A8A7

A5703021)11062A13210220图二20S6A15500420A1220A14S4S3S2

6901200720A1620211002012195306

1150

600450

23

104

A1

2.2.2运费矩阵的计算模型

问题分析

我们只考虑本题第一问的求解。首先，所有钢管必须要运到天然气

主管道铺设线路上的结点，然后才能向左或向右铺设。因此，1必须求出从每个钢管厂12

A,A,,A1妨称为运费矩阵）及其对应的运输方式和路线。因为题目中没有给出装卸成本，我们简单地假设总是采用最经济的运输方式，虽然这个假设在实际中可能不太接近事实。也就是说，在运输过程中需要多次装卸也是允许的（如铁路转公路，再转铁路，等等）。自然的想法是运输路线应该是走最短路径，但由于有两种运输和计价方式（铁路和公路），公路运输费用为1单位钢管每公里0.1万元（不足整公里部分按整公里计算），运费是路程的线性函数；然而，铁路运输要通过运输里程查表得到，是一个阶梯函数，这两种运输计价方式混合在一起，使得我们不能直接在整个铁路、公路混合的运输网络上计算最短路径作为运输路线，但可以分别在铁路、公路网上计算最短路径，然后换

dc10,

(k1)(k)ijij(kdc10,

(k1)(k)ijij(k)时的最短距离；自然u就是i，j之间的最短1d1iiu(1)

minu,uikkj可以看成任意两个接点i，j之间(n1)22，然后按照这个最w,ij,ijij(k)(k)uj,k1,2,?,i,n.求一次最短路问题，就可以把它们统一成一个标准的运费矩阵。

铁路子网络

假设铁路运输路线应该是走最短路径，而且采用连续路径计价方式

一定优于分段计价方式（其实题中数据并不符合这一假定，例如题中650km的运价计为44万元，而分成300km和350km两段计价只需要43万元，这种情况应该是不太符合实际的，可能是命题时选择数据的疏忽，我们不过多地考虑这种情况）。这时我们可以把铁路运输子网独立出来，

在这个网络上计算任意两个节点i，j之间的最短路ij短路长度查铁路运价表得到最小运费。ij在无向网络上求任意两点之间最短路的算法很多，尤其对本题这种弧上的权（距离）全为正数的情况，存在相对比较简单的算法。例如，求任意两点之间最段路径的Floyd-Warshall算法是

u(1)

u

这实际上是一种标号算法，其中n是网络中的节点数（节点编号为1，

2，„„，n）；w是给定的网络上相邻接点i，j之间的直接距离（i，jij

不相邻时取w充分大就可以了）；uijij距离的中间迭代值（或称为临时编号），既成节点i到j但不允许经过其

他节点k，k+1，„„，nij距离（或称为永久标号），即d。ij公路子网络

类似地，可以假设公路运输路线应该是走最短路径，把公路运输子网独立出来，在这个网络上计算任意两个节点i，j之间的最短路长度，

然后按照这个最短路长度乘以0.1得到最小费用c（因为公路运输费ijij用为1单位钢管每公里0.1万元）。

c1c1,S,7A,A,,A

cpsc1c1,S,7A,A,,A

cpsfxy(pc)x

y导致的运输费：这部分费用是以0.1为首项、0.1为公差的等差15jjjjj122S1y(1z)z（记为i=1，2，„„，7）到每个y

在次基础上，就可以将以上两个子网络（需要分别看成两个完全子

图）组合成一个网络每条弧上相应的运费或c为权（如果某条弧（i，ijijj上既有铁路运费，又有公路运费c，只需要取其中较小的一个即可）。ijij此时，再计算从每个钢管厂S12节点（记为j=1，2，„„，15）的最短路，得到的就是每个1215单位钢管的最小运费。ij算法仍然可以采用Floyd-Warshall算法。

2.2.3运输量的计算模型及求解

记第i个钢厂的采购费用为，最大供应量为，最小供应量为500，ii从第i个钢厂到铺设节点j的运输费用为c；用b表示管道第j段需要ijj铺设的钢管量。这些已经是已知的（或已经求得的）。

决策变量：用表示钢厂i是否使用；是从钢厂i运到节点j的iij钢管量；是从节点j向左铺设的钢管量；z是从节点j向右铺设的钢jj管量。目标函数：目标函数应该包括两部分费用。（1）从第i个钢厂采购钢管并将钢管运到铺设节点j的运费费用，对所有i，j求和，即

iijiji,j（2）从铺设节点j向左铺设的钢管量和从节点j向右铺设的钢管量j

zj级数，所以对所有j求和，可得这部分费用为

0.12

约束条件是比较明显的，所以下面直接给出本题第一问求运输量的

2

7AAyyxyz问题描述0.115

xyz,j1,2,?15,(2)152

7AAyyxyz问题描述0.115

xyz,j1,2,?15,(2)15;yiijijjjjji,jj1s.t.

500fxsf,i1,2,?7,(1)iijiij1模型如下：ijjji=1yzb,j1,2,?14,(3)j+1jjyz0,(4)115f0,1,i1,2,?7.(5)i约束（1）表示每个钢厂的生产总量限制（不低于500t，也不超过总能力限制），约束（2）表示运输到每个辅助节点的钢管数量刚好等于从这个节点向左和向右铺设的钢管的钢管数量；约束（3）表示每一段管道上铺设的钢管数量正好等于前一个节点向右、后一个节点向左铺设的

钢管数量之和；约束（4）是对端点和的自然限制；约束（5）表示115是否使用某个钢厂的0-1限制。这是一个二次规划模型。有些读者可能已经注意到，由于题目中说

明运输不足整数公里需要按照整数公里计算，所以上面费用中的第二项

只有当，z为整数时才是精确成立，否则只能看成是一种近似。不过，jj由于题目中所给的距离都是整数，所以虽然上面的模型中最优解未必一

定是整数，但这个问题必然存在，z为整数的解（自然，此时也是jjij整数），因此我们只需要在模型中再加上（或）为整数的限制就完全jj没有问题了。这样的性质通常为“占优”性质，在优化建模中利用这种性质有时是非常方便的。

2.3露天矿生产的车辆安排

2.3.1

2003高教社杯全国大学生数学建模竞赛题目B题露天矿生产的车辆安排钢铁工业是国家工业的基础之一，铁矿是钢铁工业的主要原料基地。许多现代化铁矿是露天开采的，它的生产主要是由电动铲车（以下简称电铲）装车、电动轮自卸卡车（以下简称卡车）运输来完成。提高这些大型设备的利用率是增加露天矿经济效益的首要任务。

露天矿里有若干个爆破生成的石料堆，每堆称为一个铲位，每个铲位已预先根据铁含量将石料分成矿石和岩石。一般来说，平均铁含量不低于25%的为矿石，否则为岩石。每个铲位的矿石、岩石数量，以及矿石的平均铁含量（称为品位）都是已知的。每个铲位至多能安置一台电铲，电铲的平均装车时间为5分钟。卸货地点（以下简称卸点）有卸矿石的矿石漏、2个铁路倒装场（以下简称倒装场）和卸岩石的岩石漏、岩场等，每个卸点都有各自的产量要求。从保护国家资源的角度及矿山的经济效益考虑，应该尽量把矿石按矿石卸点需要的铁含量（假设要求都为29.5%1%，称为品位限制）搭配起来送到卸点，搭配的量在一个班次（8小时）内满足品位限制即可。从长远看，卸点可以移动，但一个班次内不变。卡车的平均卸车时间为3分钟。所用卡车载重量为154吨，平均时速28。卡车的耗油量很大，每个班次每台车消耗近1吨柴油发动机点火时需要消耗相当多的电瓶能量，故一个班次中只在开始工作时点火一次。卡车在等待时所耗费的能量也是相当可观的，原则上在安排时不应发生卡车等待的情况。电铲和卸点都不能同时为两辆及两辆以上卡车服务。卡车每次都是满载运输。每个铲位到每个卸点的道路都是专用的宽60的双向车道，不会出现堵车现象，每段道路的里程都是已知的。一个班次的生产计划应该包含以下内容：出动几台电铲，分别在哪些铲位上；出动几辆卡车，分别在哪些路线上各运输多少次（因为随机因素影响，装卸时间与运输时间都不精确，所以排时计划无效，只求出各条路线上的卡车数及安排即可）。一个合格的计划要在卡车不等待条件下满足产量和质量（品位）要求，而一个好的计划还应该考虑下面两条原则之一：1.总运量（吨公里）最小，同时出动最少的卡车，从而运输成本最小；2.利用现有车辆运输，获得最大的产量（岩石产量优先；在产量相同的情况下，取总运量最小的解）。请你就两条原则分别建立数学模型，并给出一个班次生产计划的快速算法。针对下面的实例，给出具体的生产计划、相应的总运量及岩石和矿石产量。某露天矿有铲位10个，卸点5个，现有铲车7台，卡车20辆。各卸点一个班次的产量要求：矿石漏1.2万吨、倒装场Ⅰ1.3万吨、倒装场Ⅱ1.3万吨、岩石漏1.9万吨、岩场1.3万吨。铲位和卸点位置的二维示意图如下，各铲位和各卸点之间的距离（公里）如下表：铲位1铲位2铲位3铲位4铲位5铲位6铲位7铲位8铲位9铲位10

矿石漏5.265.194.214.002.952.742.461.900.641.27倒装场1.900.991.901.131.272.251.482.043.093.51Ⅰ岩场5.895.615.614.563.513.652.462.461.060.57岩石漏0.641.761.271.832.742.604.213.725.056.10倒装场4.423.863.723.162.252.810.781.621.270.50Ⅱ各铲位矿石、岩石数量(万吨)和矿石的平均铁含量如下表：

0．951．051．001．051．101．251．051．301．351．251．251．101．351．051．150．951．051．001．051．101．251．051．301．351．251．251．101．351．051．151．351．051．151．351．2530%28%29%32%31%33%32%31%33%31%C7矿石量岩石量铁含量

2.3.2运输计划模型及求解

问题分析

于2.2节的问题类似，本题也可以看成是经典运输问题的一种变形和扩展，它与典型的运输问题明显有以下不同：（1）这是运输矿石与岩石两种物资的问题；（2）属于产量大于销量的不平衡运输问题；（3）为了完成品位约束，矿石要搭配运输；（4）产地、销地均有单位时间的流量限制；（5）运输车辆只有一种，每次都是满载运输，154t/车次；（6）铲位数多于铲车数意味着最优的选择不多于7个产地作为最后结果的产地（7）最后求出各条路线上的派出车辆数及安排每个运输问题对应着一个线形规划问题，以上不同点对它的影响不同，（1）（2）（3）（4）条可通过变量的设计、调整约束条件实现；（5）条是整数要求将使其变为整数线形规划；（6）条不容易用线形模型实现，一种简单的办法是从=120个整数规划中取最优的即得到最佳物流；10为完成（7）条，由最佳物流算出各条路线上的最少派车数（整数）再给出具体安排即完成全部计算。然而这是个实际问题，为了及时指挥生产，题中要求算法是快速算法，而整数规划的本质是NPC问题，仅就20min内计算含50个变量的整数规划来说就不一定办得到。从另一个角度看，这是一个二层规划问题，第二层规划是组合优化，如果求最优解计算量较大，现成的各种算法都无能为力。于是问题变为找一个寻求近优解的近似解法，例如可用启发式方法求解。

调用120次整数规划可用三种方法避免：1先不考虑电铲数量约束运行整数规划，再对解中运量最少的几个铲位进行筛选；2在整数线形规划的铲车约束中调用sign函数来实现；3增加10个1-0个变量了来标志各个铲位是否有产量。从每个运输问题都有目标函数的角度看，这又是一个多目标规划问题，第一问的主要目标有：1重载路程最小；2总路程最小；3出动卡车数最少。仔细分析可得：1和2在第一层，3在第二层；1与2基本等价，于是只用1于第一层，对其结果在第二层中派最少的卡车，实现全局目标生产成本最小。第二问的主要目标有：4岩石产量最大；5矿石产量最大；6运量最小。三者之间的关系根据应该理解为字典序。合理的假设主要有：（1）卡车在一个班次中不应发生等待或熄火再启动的情况（2）在铲位或卸点处由两条路线以上造成的冲突问题面前，我们认为只要平均时间能完成任务，就认为不冲突。不进行具体排时方面的讨论（3）空载与重载的速度都是28km/h，耗油相差却很大（4）卡车可提前退出系统，等等

模型建立

我们将只考虑本题第一问（求运输量）的模型。首先定义如下数学符号：Xij为i号铲位到j号卸点路线上运行一个周期平均所需要时间（min）；Aij为从i号铲位到j号卸点最多能同时运行的卡车数（辆）；Bij为从i号铲位到j号卸点路线上一辆车最多可以运行的次数（次）；Pi为i号铲位的矿石铁含量（%）；P=（3028293231332323331）；qj为i号卸点任务需要（t）；q=（1.2，1.3，1.3，1.9，1.3）*10000。CKi为i号铲位的铁矿石储量（万t）；CYi为i号铲位的岩石储量（万t）；fi为描述第i号铲位是否使用的0-1开关变量，取1为使用，取0为关闭。目标函数与约束条件的分析；（1）目标函数取为重载运输时的运量（.t.km）最小（因Xij代表整数车次数，乘154后等于运量；再乘以运输距离等于吨公里）。（2）道路能力约束；由于一个电铲（卸点）不能同时为两辆卡车服务，所以一条路线上最多能同时运行的卡车数是有限的，卡车在i号铲位到j号卸点路线上运行一个周期平均所需的时间为T=（i到j距离*2/平均速度）+3+5（min）由于装车时间5min大于卸车时间3min，所以可分析出这条路线上在卡车不等待条件下最多能同时运行的卡车数为：Aij=[Tij/5]。同样可分析每辆卡车一个班次中在这条路线上最多可以运行的次数为Bij=[8*60-（Aij-1）*5]/Tij，其中（Aij-1）*5是开始装车时最后一辆车的延时时间。则一个班次中这条固定路线上最多可能运行的总车次大约为：Lij=Aij*Bij。（3）电铲能力约束：还是因为一台电铲不能同时为两辆卡车服务，所以一台电铲早一个班次中的最大可能产量为：8*60/5

1545

10

10

1545

10

10

10

10

10

xij

ij

T105xf*8*60/5,i1,2,…，10，

x8*20,j1,2,

xq/154,j1,2,，

x*(p30.5)0,j1,2,

x*(p28.5)0,

f7,

20,

ij5（4）卸点能力约束：卸点的最大吞吐量为每小时60/3=20车次。于是一个卸点在一个班次中的最大可能产量为：8*20车（5）铲位储量约束：铲位的矿石和岩石产量都不能超过相应储藏量（6）产量任务约束：各卸点的产量大于等于该卸点的任务要求（7）铁含量约束：各矿石卸点的平均品位要求都在指定的范围内（8）电铲数量约束：电铲数量约束无法用普通不等式表达，但算法运行时间增加了，也可以用其他方法解决该问题（9）卡车数量约束：卡车总数不超过20辆（10）整数约束：当把问题作为整数规划模型时，车流量Xij为非负整数这样，我们可以得到的各种模型之一为：

minxc,ijiji1j1s.t.xAB,i1,2,…，10，j1,2,„，5，ijijijijij1iji1xxxck*10000/154,i1,2,…，10，i1i2i5ixxcy*10000/154,i3i4iijji1ijii1ijii1ii1

i,f0,1;x为整数，i1,2,„，10，j1,2,„，5.iij其中，Ti到j距离*2/平均速度3（5min），ij

A,B8*60(A1)*5/T近似）ij=ijijij

空洞探测PQPQt

i

PQtQjPiij0.05920.41310.44530.4529空洞探测PQPQt

i

PQtQjPiij0.05920.41310.44530.45290.31770.3397jijCBQ20.08950.44130.05980.40400.22630.23640.3566Q30.19960.43180.41530.07380.19170.30640.1954Q40.20320.47700.41560.17890.08390.22170.0760Q50.41810.52420.35630.07400.17680.09390.0688Q60.49230.38050.19190.21220.18100.10310.1042Q70.5646

2.4.1问题描述

山体、隧洞、坝体等的某些内部结构可用弹性波测量来确定。一个简化问题可描述为，一块均匀介质构成的矩形平板内有一些充满空气的空洞，在平板的两个邻边分别等距地设置若干波源，在它们的对边对等地安放同样多的接收器，记录弹性波由每个波源到达对边上每个接收器的时间，根据弹性波在介质中和在空气中不同的传播速度，来确定板内空洞的位置。现考察如下的具体问题：一块240（米）×240（米）的平板（如图），在AB边等距地设置7个波源(i=1,„,7)，CD边对等地安放7个接收器(j=1,„,7)，记录由发出的弹性波到达i的时间(秒);在AD边等距地设置j7个波源R(ii=1,„,7)，BC边对等地安放j7ij个接收器S(j=1,„,7)，记录由R发出的弹性波到达iS的时间τ(秒)。已知弹性波在介质和空气中的传j播速度分别为2880（米j/秒）和ij320（米/秒），且弹性波沿板边缘的传播速度与在介质中的传播速度相同。1）确定该平板内空洞的位置。2）只根据由发出的弹性波到达的时间(i,j=1,„,7)，能确定空洞的位置吗；讨论在同样能够确定空洞位置的前提下，减少波源和i

D

Ri

Sj

A接受器的方法。

tQ1P0.0611P10.0989P20.3052P30.3221P40.3490P50.3807P670.4311

ij0.07530.34560.36550.31650.27490.443440m40m的正方形平板被均匀地分成个小正方形单元，每个，k,l6xxS0.064510.07000.32050.32890.24090.3891ij0.07530.34560.36550.31650.27490.443440m40m的正方形平板被均匀地分成个小正方形单元，每个，k,l6xxS0.064510.07000.32050.32890.24090.38910.4919例S0.060220.28520.09740.42470.32140.58950.3904如S0.081330.43410.40930.10070.32560.30160.0786单S0.351640.34910.42400.32490.09040.20580.0709元S0.386750.48000.45400.21340.18740.08410.0914（S0.431460.49800.31120.10170.21300.07060.0583k,lS0.57217）表示由点（40RR1R2R3R4R5R67

2．4.2优化模型及求解

问题分析

这道题目有多种解法。测量问题一般是随机性问题，我们把这个问题看成是确定性问题来解。虽然平板上空洞的大小和形状可以是任意的，我们这里把它简化成的正方形（可以成为单元），即240m240m6636小正方形或者全是空洞，此外，在以下假设下考虑这个问题的解法：（1）观测数据有测量误差，观测数据除测量误差外是可靠的。（2）波在传播过程中沿直线单向传播，且不考虑波的反射，折射以及干涉等现象。（3）空气和介质都是均匀的。（4）“弹性波”在传播过程中没有能量损失，其波速仅与介质有关，且在同一均匀介质中波速不变。（5）假设平板可划分为网络，空洞定位于每个网格单元内，空洞大小大致相同。（6）题中已经假设弹性波沿板边缘的传播速度与在介质中的传播速度相同，此外，网格化以后，如果某条波线位于平板中两个单元的边缘（交线）上，我们假设这条波线的传播是沿着每个单元各走一半的路程（或者换句话说，如果这两个单元中一个是空洞，另一个是介质，则波线沿这条交线的传播速度是在空洞和介质中传播速度的平均值）。做这样的网格化（离散化）以后，设一边上的波源为m个（本题m=7个），另一边上的波源为n个（本题n=7）,则得到(m-1)(n-1)个单元。

建立平面直角坐标系，设A点为原点，AB为x轴，AD为y轴，则C点坐标为（240，240），假设平板上每个单元是从左到右。从下到上编号的（k-1）,40(l-1)）,(40k,40(l-1),(40(k-1),40l),(40k,40l)围成的小正方形单元，其中k,l的取值范围是1。每个单元对应一个决策

变量，表示该单元是否为空洞（1表示是，0表示否）。我们的任务是kl确定这36个0-1变量的取值。为此，首先需要计算每条波线与每个单kl元的交线的长度，这虽然可以使用任何软件编程进行计算，但我们下面还是说明如何用LINGO进行计算。

Rk,l7QPQjQPRk,l7QPQjQPQ）Q240(ki)/(ji),其中l16(ki)/(ji)l16(ki)/(ji)l40,如果kij1或k1ij7;20,y=6(x-40(i-1)).1如果kij或k1ij,且2k5l6

可以先考虑波源P与接受器Q决定的波线与每个单元（k,l）的交线，ij记其长度为b（对与S可以类似的考虑，或直接根据对称性得到结果，ijklij我们将在后面再讨论）。假设P与Q都是从左到右编号的，其中i,j的取ij值范围是1，如P位于原点（0，0），Q位于（40，240），等等。12如果i=j，则波线P与y轴平行，此时只要不是平板上最左边和最ij右边的波线，波线都会位于两个单元的边缘（交线）上，根据上面ij的假设6，容易得出：

bijkl（43）

下面来考虑i的情况。

因为波源P的坐标是（40(i-1),0）,的坐标是（40（j-1）,240）,ij所以容易得到的直线方程为ij（j-i(44)这条直线P与每个单元（k,l）的边缘最多只能有两个交点，但交ij点有可能位于单元的不同的边缘位置（每个单元有上，下，左，右四个边缘位置）。虽然对于像本题这样规模不太大的问题，可以通过枚举法确定所有交点，但我们下面还是介绍能够用于更大规模问题的一般解题方法。对于单元（k,l）的左边缘，其对应的直线方程为

x=40(k-1).(45)将（45）代入（44），可以得到波线与单元边缘对应交点的y坐标为

y1ijkl（46）式（46）条件l是为了保证这个交点是有效的，

即这个交点确实位于这个单元（k,l）所在的范围内。

(120i20(ij)(l1))(120i20(ij)(l1))/3k1)(120i20(ij)(l1))/340k06(ik)(ij)(l1)6，其中(120i20(ij)l)/340(k1)(120i20(ij)l1)/340k6(ik)(ij)l640l,其中06(ik)(ij)l6Q

PQj06(ik)(ij)(l1)6。

x40k(47)将（47）代入（44），可以得到对应交点的y坐标为

y2240(k1i)/(ji),其中l16(k1i)/(ji)lijkl（48）对应单元（k,l）的下边缘，其对应的直线方程为

y40(l1)（49）将（49）代入（44），可以得到对应交点的x的坐标为

x3ijkl（50）但（50）只有在40(时才是有效的，

这个条件化简后就是，于是可以得到对应交点

的y坐标为

y340(l1)ijkl（51）同理，对于单元（k,l）的上边缘，其对应的子痫方程为

y40l（52）将（52）代入（42），可以得到对应交点的x的坐标为

x4ijkl（53）但是（53）只有在时才是有效

的，这个条件化简后就是0，于是可以得到对应的交

点的y坐标为

y4ijkl（54）至此，式（46），（48），（51），（54）给出了直线P与单元（k,l）ij的边缘所有可能的交点的y坐标及其存在交点的条件，但事实上，这些条件