运筹学模型与数学建模竞赛

一、引言

一般来说，大学生数学建模竞赛所涉及到的运筹学模型包括数学规划（线性规划和非线性规划），网络优化（含网络计划技术），排队模型，动态规划等，请看下表

下面重点介绍运筹学模型的数学规划。二、数学规划的一般形式

))(m ax ()(m in x f or x f

??≤≤=≤==ub x lb m j x g l i x h t s j i ,,2,1,

0)(,,2,1,0)(.

. 线性规划：整数规划：非线性规划：

三、数学规划问题举例

1 下料问题

现要用100×50厘米的板料裁剪出规格分别为40×40 厘米与50×20厘米的零件，前者需要25件，后者需要30件。问如何裁剪，才能最省料？

解：先设计几个裁剪方案

记 A---------40×40

注：还有别的方案吗？

显然，若只用其中一个方案，都不是最省料的方法。最佳方法应是三个方案的优化组合。设方案i 使用原材料x i 件(i =1,2,3)。共用原材料f 件。则根据题意，可用如下数学式子表示：

???

??=≥≥++≥+++=)3,2,1(0305325

.min 321213

21j x x x x x x t s x x x f j

，整数这是一个整数线性规划模型。

2 运输问题

现要从两个仓库（发点）运送库存原棉来满足三个纺织厂（收点）的需要，数据如下表，试问在保证各纺织厂的需求都得到满足的条件下应采取哪个运输方案，才能使总运费达到最小？（运价（元/吨）如下表）

方案1

方案2

方案3

解：题意即要确定从i 号仓库运到j 号工厂的原棉数量。故设ij x 表示从i 号仓运到j 号工厂的原棉数量（吨）f 表示总运费.则运输模型为：

??????

?????

??==≥?

??=+=+=+???

≤++≤+++++++=运输量非负约束；需求量约束运出量受存量约束),,j ,i (x x x x x x x x x x x x x .t .s x x x x x x f min ij 32121025154030504223223

1322122111232221131211232221131211 一般地，对于有m 个发点和n 个收点的运输模型为

?????

??????==≥===≤=∑∑∑∑====),...2,1;,...2,1(0)...2,1(),...3,2,1(..min 1

n j m i x n j b x m i a x t s x c f ij m

i j ij n

j i ij m i n

j ij

ij 其中a i 为i 号发点的运出量，b j 为j 号收点的需求量，c ij 为从i 号发点到j 号收点的单位运

价。特别当

∑∑===m i n

j j

i b

a 1

时，存货必须全部运走，故上述约束条件中的

∑=≤n

j i ij

a x

可改为等式：

),...2,1(1

m i a x

i n

j ij

==∑=

3 选址问题

某地区有m 座煤矿，i # 矿每年产量为a i 吨，现有火力发电厂一个，每年需用煤b 0吨，每年运行的固定费用（包括折旧费，但不包括煤的运费）为h 0元。现规划新建一个发电厂，m 座煤矿每年开采的原煤将全部供给这两个电厂发电用。现有n 个备选的厂址。若在j #备选厂址建电厂，每年运行的固定费用为h j 元，每吨原煤从i # 矿运送到j #备选厂址的运费为c ij 元(i =1,2,…m , j =1,2…n )。每吨原煤从i # 矿运送到原有电厂的运费为c i0 (i =1,2,…m )。试问：

[1] 应把新电厂厂址选在何处？

[2] m 座煤矿开采的原煤应如何分配给两个电厂？

才能使每年的总费用（电厂运行的固定费用与原煤运费之和）为最小？

模型的建立

(1) 变量的设置为了解决问题[1]，我们使用0-1变量

n j j y j ,2,10

#1=??

?=否则

备选厂址选中

为了解决问题[2]，设从i #煤矿运到j #备选的厂址的运量为x ij 吨（i=1,2,…m ，j=0,1,2,…,n ）

总运费：

ij m i n

j ij

x c

∑∑==1（对不被选中的备选厂址运费x ij,将由约束条件限制为0）.

固定费用 h 0+

∑=n

j j j

y h

每年总费用 z =

h y h x c

j j j m

i n

j ij ij

++∑∑∑===

(3)约束条件的表达（i ）煤矿产量约束

m ,,i a x

j ij

210

==∑=

（ii ）旧电厂用煤量约束

b x

i i =∑=

（iii ）新电厂用煤量约束记 01

b a

b m

i i

∑=，当j #

备选厂址被选中时∑==m

i ij b x 1

，当

j #

备选厂址没被选中时

∑==m

i ij

0，综合表达为n j y b x j

i ij ,...2,11

==∑=

（iv ）选址约束由于只选一个厂址，所以

=∑=n

j j

(v)非负及整数约束

j y n

j m

i x j ij ,2,11

0,2,1,0,2,10====≥或

综合得数学规划模型：

00011

min ,1,...,,1,...,..10,1,...,;0,...,0,1;1,...,m n n

ij ij j j i j j n

ij i j m

i i m

ij j i m i i n

j j ij

j z c x h y h x a i m x b x by j n s t b a b y x i m j n y j n =========++?==???=???==????

=-???=???≥==??==????

∑∑∑∑∑∑∑∑ 4布点问题

某市有6个区，每个区都可建消防站，为了节省开支，市政府希望设置的消防站最少，但必须保证在该市任何地区发生火警时，消防车能在15分钟内赶到现场。假定各区的消防站要建的话,就建在区的中心，根据实地测量,各区之间消防车行使的最长时间如下表:(单位:分钟)

请你为该市制定一个设置消防站的最节省的计划。建模并求解。

解：本题实际上是要确定各个区是否要建立消防站，使其既满足要求，又最节省。这自然可引入0-1变量，故设

）

（区建消防站时

当不在第区建消防站时

当在第621 0 1,,,j j ,j ,x j =???= 目标是∑==

j j

f 最少。以下考虑约束条件。

若1区发生火警，按照“消防车要在15分钟内赶到现场”的要求，则l 区和2区至少应有一个消防站，即121≥+x x 。同理得：

1 1 1 1 165265454343621≥++≥++≥++≥+≥++x x x ,x x x ,x x x ,x x ,x x x

从而得模型为：

621101*********

5454343621216

?????

???

??==≥++≥++≥++≥+≥++≥+=∑=)

,,,j (,,x x x x x x x x x x x x x x x x x .t .s x f min j j j

再仔细观察知，若满足第一、三个约束条件，则必然满足第二、四个约束条件，故后者是多余的，可省略。从而化简得：

621101

5265443216

?????

??==≥++≥++≥+≥+=∑=)

,,,j (,,x x x x x x x x x x x .t .s x f min j j j

此模型当然可用软件求解，但由于比较简单,故可直接试算。若要求只有一个1=j x ，则显然不可行；若要求只有两个1=j x ，则有唯一可行解0 1653142======x x x x ,x x ，故这就是最优解，即只需在2区和4区建立消防站。

附程序：

c=[1 1 1 1 1 1]

a=-[1 1 0 0 0 0;0 0 1 1 0 0;0 0 0 1 1 1;0 1 0 0 1 1]

b=-[1 1 1 1]

v=zeros(1,6)

u=[1 1 1 1 1 1]

[x fval]=linprog(c,a,b,[],[],v,u)

Optimization terminated.

x =

0.0000

1.0000

0.0000

1.0000

0.0000

fval =

2.0000

5 配套问题

设有n个车间要生产m种产品，第j车间每天生产第i种产品至多a ij件（即全天只安排生产产品i而不安排生产其他产品时的最大产量），假设这m种产品每种一件配成一套，问如何安排生产任务才能使生产出的成套产品最多？（i=1,2,...,m；j=1,2,...,n）

（1）建模方法（一）

设x ij ——车间j 安排用于生产产品i 的数量，Z ——每天生产的成套产品数目，

原问题可以转化为以下数学模型：

max min n

ij i m

j f z x ≤≤===∑

..0(1,

,;1,,)

ij ij

ij f z s t x a x i m j n ?≥?≤??

≥ ==?

模型改进为：

max f z =

1(1,,)..0(1,,;1,,)

ij j ij ij

ij x z i m s t x a x i m j n =?≥ =???≤??

≥ ==???

∑

模型问题：没有将一天的生产时间约束考虑到！

2、建模方法（二）

设 x ij ——车间j 安排用于生产产品i 的时间（占全天的比例）

Z ——每天生产的成套产品数目

则a ij x ij ——车间j 每天生产产品i 的数目。例如：车间2每天至多生产某产品6件，若安排1/3天时间去生产，则可产出2件），而

a ij

j ij

∑=1

——每天全厂产出产品i 的总量。

于是则有模型

max Z ( 11min n

ij i m

j z x ≤≤==∑)

(*) ??????

?????≥==≥=≤=≥∑∑==整数

02121 0211211

Z )n ,...,,j ;m ,...,,i (x )n ,...,,j ()m ,...,,i (Z .t .s ij m i ij n

j ij ij x x a

其中常数1表示1天。

注：(1)此模型着重考虑安排生产的时间；(2)从实际情况考虑，安排生产的时间必须是每件产品耗用生产时间的整数倍才合适。例如，每3分钟生产一件产品，安排5分钟，也只能生产1件，不能认为生产了1.67件。模型（*）没有考虑到这些因素，故是不合适的。

（2）建模方法（二）——改进（*）

显然，

a 1

——车间j 生产每件产品i 的耗用时间（天）。从以上分析应有 ?

??=a x ij ij 1? (?是非负整数)

从而令 y ij =a ij x ij , y ij 是非负整数,表示车间j 每天生产产品i 的数目，将它代入（*）后得

max Z

(**) ?????

???

≥==≥=≤???

? ??=≥∑∑==整数整数02121021112111

Z )n ,...,,j ;m ,...,,i ()n ,...,,j ()m ,...,,i (Z .t .s y y a y ij m i ij ij

j ij 这是一个整数线性规划模型。

注：此模型着重考虑安排生产产品的数目。

四、数学规划在数学建模中的应用举例

1．钻井布局

勘探部门在某地区找矿。初步勘探时期已零散地在若干位置上钻井，取得了地质资料。

进入系统勘探时期后，要在一个区域内按纵横等距的网格点来布置井位，进行“撒网式”全面钻探。由于钻一口井的费用很高，如果新设计的井位与原有井位重合（或相当接近），便可利用旧井的地质资料，不必打这口新井。因此，应该尽量利用旧井，少打新井，以节约钻探费用。比如钻一口新井的费用为500万元，利用旧井资料的费用为10万元，则利用一口旧井就节约费用490万元.

设平面上有n 个点P i ，其坐标为(a i ,b i ),i =1,2,…,n ，表示已有的n 个井位。新布置的井位是一个正方形网格N 的所有结点（所谓“正方形网格”是指每个格子都是正方形的网格；结点是指纵线和横线的交叉点）。假定每个格子的边长（井位的纵横间距）都是1单位（比如100米）。整个网格是可以在平面上任意移动的。若一个已知点P i 与某个网格结点X i 的距离不超过给定误差ε（=0.05单位），则认为P i 处的旧井资料可以利用，不必在结点X i 处打新井。

为进行辅助决策，勘探部门要求我们研究如下问题：

1)假定网格的横向和纵向是固定的（比如东西向和南北向），并规定两点间的距离为其横向距离（横坐标之差绝对值）及纵向距离（纵坐标之差绝对值）的最大值。在平面上平行

移动网格N，使可利用的旧井数尽可能大。试提供数值计算方法，并对下面的数值例子用计算机进行计算。

2)在欧氏距离的误差意义下，考虑网格的横向和纵向不固定（可以旋转）的情形，给出算法及计算结果。

3)如果有n口旧井，给出判定这些井均可利用的条件和算法（你可以任意选定一种距离）。

数值例子n=12个点的坐标如下表所示：

①问题重述（略）

②问题分析

布局问题是在一定约束条件下的最优化问题，勘探部门要求尽可能多地利用旧井，因此我们必须围绕这个问题来移动网格。

问题（1）：网格的横向和纵向是固定的，网格只能平移。如果考虑网格上所有的结点通过平移后与旧井位的距离，由于网格结点数较多，运算量较大，虽然可以解决问题，但并非一个好的解法。因此我们从旧井位考虑，通过四舍五入法找到与其相邻的网格结点（相互关系如图所示），然后我们通过两个控制变量（角度与距离）对这些结点进行平移，使网格结点尽可能多的与旧井位重合。

问题（2）：网格不固定移动，需要我们考虑平移和旋转，对于旋转我们同样从旧井位出发，通过旧井位与网格结点的相对关系，不难得出：旧井位的旋转是网格旋转的逆过程，因此我们利用旧井位旋转来替代网格的旋转问题，然后，根据问题（1）的方法，再解决网格平移问题，这样就将网格的移动问题简化成为旧井位的旋转问题和网格的平移问题。

问题（3）：以上两个问题都是在正方形网格情况下进行的求解，在此问中我们不妨假设网格是长方形的，通过调整长和宽的比例关系使n个旧井位都可以被利用。因此问题（3），就是在问题（2）的基础上，使n个旧井位都可以被利用的长方形的长、宽比例问题。具体实现方法如下：把目标函数值定为12，通过对问题（2）的反向求解来确定长方形的长、宽比例。

③模型假设和符号说明

模型假设：

1．若一个旧井位与某个网格结点不超过0.05单位，则视为重合；

2．网格是足够大且平铺的；

3．模型采用直角坐标系，在网格未移动时，网格的横向和纵向就是直角坐标系的两坐标轴方向。

4．逆时针旋转为正，顺时针方向为负；

5．长方形网格的长是1单位。

符号说明：

i...i旧井位,i=1,2, (12)

Pi…第i个旧井位;

ai,bi…第i个旧井位的横坐标、纵坐标;

Ai,Bi…第i个旧井位旋转的横坐标、纵坐标;

W i…与第i个旧井位相邻的网格结点;

xi,yi…与第i 个旧井位相邻的网格结点的横坐标、纵坐标;

xli,yli…与第i 个旧井位相邻的网格结点移动后的横坐标、纵坐标; di…第i 个旧井位与其相邻的网格结点距离; s…网格平移长度; an…网格平移角度;

bn…旧井位点旋转角度;

pi…第i 个旧井位到坐标原点的距离; qi…第i 个旧井位与坐标原点连线的倾角; T …可被利用的旧井位数;

ti…第i 个旧井位数是否被利用的函数(1表示可被利用;0表示不可被利用); c…长方形的宽与长的比例;

um, tm…分别为宽与长之比的上、下线. ④模型的建立

我们的目标是使旧井位利用数最多的T, 它由ti(d)决定，即T=

∑=12

)(i d ti 。而由要求知

ti(d)=otherwise

i d 05.0)(00

1≤≤??

? (1)

d(i)的表达式在下面给出。

由于问题（1）是问题（2）的特殊情况，所以这里只给出问题（2）的模型。我们先将旧井位进行旋转。

第i 个旧井位到坐标原点的距离为

p(i)=2

)()(i b i a +. (2)

第i 个旧井位与坐标原点连线的倾角为

q(i)=arctan(b(i)/a(i)). (3)

所以第i 个旧井位旋转的横坐标、纵坐标分别为

A(i)=p(i)cos[q(i)+bn(i)], B(i)=p(i)sin[q(i)+bn(i)] (4)

(bn 旧井位点旋转角度)。然后用四舍五入法求旧井位的相邻网格结点坐标：

x(i)=[A(i)+0.5], y(i)=[B(i)+0.5] (5)

由此得出网格结点平移后的坐标为

xl(i)=x(i)+s*cos(an), yl(i)=y(i)+s*sin(an) (6)

所以，网格结点与旧井位的误差距离为

d(i)=2

)]()([)]()([i B i yl i A i xl -+- (7)

综上所述，本问题的数学模型可表示为：

目标函数： maxT= ∑=12

)(i d ti

约束条件是满足(1)-(7)式。

注：以上模型说明：当bn=0, d(i)=max{|xl(i)-A(i)|, |yl(i)-B(i)|}时，目标函数的解即

为问题（1）的解。

⑤模型的求解

我们的目标是使目标函数最大化，这是一个有约束的优化问题，我们采用穷举法来求它的极值。所谓穷举法就是逐点确定寻优方向，然后利用固定的步长，进行搜索的方法。为使目标函数值最大，我们列出主要步骤如下：

1．定一个能包含所有解的初始范围与固定的搜索步；

2．据运行时间，再对搜索步长用穷举法进行精度调整。

根据这种方法我们得到了问题答案：

问题（1）：最多可利用的旧井位数t=4, 网格的平移方向an0=5.400000(弧度)，网格的平移距离s0=0.583000(单位)；

问题（2）：最多可利用的旧井位数t=6, 网格的平移方向an0=2.760000(弧度)，网格的平移距离s0=0.230000(单位), 网格旋转角度：bn0=-0.780000(弧度)（旧井点的相对旋转角度为：5.500000弧度）。

⑥模型的改进

降落伞的选择

一、问题的重述：

为向灾区空投一批救灾物资共2000kg，需选购一些降落伞。已知空投高度为500m，要求降落伞落地时的速度不能超过20m/s，降落伞的伞面是半径为r的半球面，用每根长L 共16根绳索连接的重m位于球心正下方球面处，如下图：

每个降落伞的价格由三部分组成，伞面费用C1由伞的半径r决定，见下表；绳索费用

面积的乘积成正比。为了确定阻力系数，用半径r=3m，载重m=300kg的降落伞从500m高

试确定降落伞的选购方案，即共需多少个伞，每个伞的半径多大（在给定半径的伞中选），在满足空投要求的条件下，使费用最低。

二、问题的分析：

根据题意，每种伞的价格是确定的。要确定降落伞的选购方案，即共需多少个伞，每个伞的半径多大（在给定半径的伞中选），在满足空投要求的条件下，使费用最低。首先，必须知道每种伞在满足空投要求的条件下最大的载重量M(r)，然后就是一个线性整数规划了。欲得到M(r)，必须先求出空气阻力系数k，然后根据运动方程得出M(r)。最后运用分支定界法求解线性整数规划，得出问题要求的解。

三、基本假设：

1．救灾物资2000kg可任意分割；

2．降落伞落地时的速度不能超过20m/s；

3．降落伞与绳索的质量可以忽略；

4．降落伞下落过程中，只受到重力和空气阻力的作用； 5．空气阻力的阻力系数k 是定值，与其他因素无关。

四、符号说明：

M (r )------ 半径为r 的伞在满足空投要求的条件下最大的载重量； t ------ 降落伞从开始下降开始记时的时间； k ------空气阻力系数

H (t )------降落伞从降落位置到t 时刻所下降的距离； m ------降落伞负重重量； g------重力加速度； s------降落伞伞面面积；

n r ------选购的半径为r 的降落伞的个数

五、模型的建立与求解

1、计算每种伞的单价如下：单位为元，C 2=4162??r

2、定空气阻力系数k 。

对给定的r=3m ，m=300kg ，取g=9.8m/s 2，s=2πr 2，有数据

处理得：

作出H(m)~t(s)的关系图1:

图1

从图1可看出：一方面，H(m)~t(s)在后阶段基本是线性关系，即降落伞作匀速运动。从

mg =kVs

有又由于若在500米的物体做匀速运动，则50017/30

s v m s t =

=，代回mg =kVs 中，估算出 k =2.9。

另一方面，根据题意，物体在降落过程中一直做匀速直线运动是不可能的。

题中注意到：降落伞在下降过程中只受到重力和空气阻力的作用，而且初速度为0，由牛顿第二定理知：

()dV t F ma m

mg kvs dt ===-，()dV t kvs

g dt m

? =-

?????=-=0

)0()

(V m kvs g dt t dV （1）解得: ks

mge

ks mg t V m

kst

--=

)( （2）

由（2）式，并代入k =2.9，r =3m ，m =300kg ，取g=9.8m/s 2，s=2πr 2，作出下图2

图2

从图2可以看出当09t ≤≤时，V (t)迅速增长，但由于负项是成负指数衰减的，所以很

快就接近极限值mg/ks 。当9t ≥时，0,kst m

kst m

mge mg ks

kse

→(),mg

v t ks

→

所以9秒以后可以看作近似的匀速运动。下面就9秒以后的数据运用最小二乘法进行线性拟合。

设H (t)=V t+b+δ，其中δ服从正态分布。Matlab 程序如下： x=[9 12 15 18 21 24 27 30];

H=[128 183 236 285 340 392 445 499]; P=polyfit(x,H,1)

从而得到p=[17.5794 –29.2976]，所以V=17.5794(m/s)，并由mg =kVs 得到k=2.9575。

3、瞬时速度与高度

设降落伞从降落位置到t 时刻所下降的距离为H (t)，则有

（3）

积分求得

2222

)(s

k g

m s k ge

m ks mgt t H m

kst

（4）

4、求半径为r 的降落伞在满足空投要求的条件下最大的载重量M (r ).

V(t)

V(t) S 0 t 0 t S=H(t 0)=0

()t v t dt ?

2222

()()kst m

kst

m mg mge V t ks ks

mgt m ge m g H t ks k s k s --

??=-????=+-?

()v t 与t 直接相关，也与m 相关。m 越大，则v 越大。由（2）式知：V(m)是关于m 的增函

数。当v 最大时，m 也达到最大M ，即为最大的载重量。

特别的，在给定从500m 的高空空投时，降落伞落地瞬间的速度在给定g, k, s 后又有等式约束H(t)=500，即

2222

500s

k g m s k ge m ks mgt m

kst

的情况下，V(m)是关于m 的增函数；反之，其反函数也是关于V 的增函数。所以要求半径为r 的降落伞在满足空投要求的条件下最大的载重量M (r )，就是要在V 取最大值时取得，即取V=20m/s 时求出指定半径r 的M(r)，于是得到如下方程组：

mge

ks mg t V m

kst

)(

2222

22)(s

k g

m s k ge

m ks mgt t H m

kst -+

=- （5）由此导出：

s k mg ksV g m t H --

?-=)/()1ln()(222 （6）如前所述，取H(t)=500，V=20m/s ，得到方程：

s k mg ksV g m --

?-=)/()1ln(500222 （7）将参数g=9.8，k=2.9575代入上式，有

s m s m ?-

??-

?-=9575.220)9575.2/()8.99575.2201ln(8.9500222 由s=2πr 2，分别代入r=2，r=2.5；r=3；r=3.5；r=4，并调用Matlab 的命令solve ，分别解

得半径为r 的降落伞在满足空投要求的条件下最大的载重量M (r )如下：

3． 4．原问题就是如下的一个线性整数规划问题：

设n r 为选购的半径为r 的降落伞的个数，则有

}15621177822596446m in{45.335.22n n n n n ++++

s.t. 2 2.53 3.541512373414646062000

,2,2.5,3,3.5,4r r n n n n n n r n N ++++≥??

=∈?

（8）

运用分支定界法可求得：6,

0345.35.22=====n n n n n

总费用为4932元，即在需要空投2000kg 的情况下，需要选购半径为3米的降落伞6把。

六、模型的检验及推广

1．在求解空气阻力系数k 时，我们分析数据得出在运动后期降落伞作近似的匀速运动，并由此为前提对数据进行拟合求出了k 。下面在问题已经基本得到解决后，运用获得数据对降

可以看出，M/s 几乎是常数，又有mg =kVs 是降落伞后期运动为匀速运动的充分必要条件，即 m/s=kV/g 为常数，而k ，g 为常数，所以降落伞在运动后期为近似匀速运动。因此，在求解空气阻力系数k 时，假设后期运动为近似匀速运动是合理的。

2．从图2可以看出降落伞的大幅度加速过程很快就结束了，对给定的伞和给定的承载质量很快就进入近似匀速运动，而且速度与空投高度基本上是无关的，所以空投高度并不十分重要，只要能保证空投位置准确，高一些投放也可以，这样可降低空投高度。

3．题目中要求的是空投2000kg 的救灾物资，若数据有变动，例如3000kg ，5000kg ，则只需将（8）式的第一个约束右端改为3000kg ，5000kg ，然后再运用分支定界法可求得结果。