运筹学模型与数学建模竞赛
1、引言
一般来说,大学生数学建模竞赛所涉及到的运筹学模型包括数学规划(线性规划和非线性规划),网络优化(含网络计划技术),排队模型,动态规划等,请看下表
注:从年起,全国大学生数学建模竞赛开始设置专供大专院校学生做的题。
下而重点介绍运筹学模型的数学规划。
二、数学规划的一般形式
nin f(x) (ornnx f(x))
/l, (x) = 0, i = 1,2,…丿
s.t.<0, ) = 12…,加
lb 线性规划: 整数规划: 非线性规划: 三、数学规划问题举例 1下料问题 现要用100X50厘米的板料裁剪出规格分别为40X40厘米与50X20厘米的零件, 前者需要25件,后者需要30件。问如何裁剪,才能最省料? 解:题意即要确立从i 号仓库运到j 号工厂的原棉数量。故设X”表示从i 号仓运到j 号工厂的原棉数量(吨)f 表示总运费?则运输模型为: min f = 2x H +X|2 +3^13 + 2x 2| + 2x 22 + 4x 23 ■ x H +x [2 +X13 S 50 x 21 + x 22 + x 23 < 30 X 11+X 2I =40' s 』:X [2+X 22 =15需求量约束 + AS j =25 列no 仃= 1,2;丿? = 123丿运输量非负约束 一般地,对于有m 个发点和门个收点的运输模型为 n 工? 5q(7 = h2,3,??m) m /=i Xq nO(j = 12??〃;J = 12??n) 其中q 为i 号发点的运出量,bj 为j 号收点的需求咼,5为从i 号发点到j 号收点的单位运 价。 m n n 特别当工% =工耳时,存货必须全部运走.故上述约朿条件中的工耳可改为等式: r-1 j-1 n 工七=£(,= 1,2,...w ) 3选址问题 某地区有m 座煤矿,尸矿每年产量为q 吨,现有火力发电厂一个,每年需用煤b 。吨, 每年运行的固左费用(包括折旧费,但不包括煤的运费)为ho 元。现规划新建一个发电厂, m 座煤矿每年开采的原煤将全部供给这两个电厂发电用。现有门个备选的厂址。若在尸备 选厂址建电厂,每年运行的固左费用为%元,每吨原煤从严矿运送到严备选厂址的运费为 5元(口j=1,2 -n )o 每吨原煤从厂矿运送到原有电厂的运费为细(i=1,2,...m )。 试问: [1] 应把新电厂厂址选在何处? [2] m 座煤矿开采的原煤应如何分配给两个电厂? 才能使每年的总费用(电厂运行的固左费用与原煤运费之和)为最小? 运岀量受存量约束 min m n f = H C u X U 模型的建立 (1)变量的设置为了解决问题[1],我们使用0?1变量 [1选电/#备选厂址.“ y; =< — , / = 1,2< -n 70否则 目标函数的表达 m n 总运费工工勺勺(对不被选中的备选厂址运费x场将由约束条件限制为0)?r-l 每年总费用SW6jXij + 22h jy J +仏z-l J—O y-1 (3)约束条件的表达 (i)煤矿产量约束 n 工? = q i = 12…加m (ii)旧电厂用煤量约束2>,。=九 (iii)新电厂用煤疑约束m m 记b =工山_叽,当严签选厂址被选中时工切=〃,当 r-l r-l tn m 严备选厂址没被选中时工切=0,综合表达为工勺=b y) j = 1,2..../? >-1 J-I (iv)选址约束由于只选一个厂址,所以丈=1 7-1 x.j > 0 i = 1,2,…m j = 0丄2, - n (V)非负及整数约束 y} =0或1 j = 12 综合得数学规划模型: min乙=为为°护j +为山儿+心 /=! >=() ./=! 艺?=曲=1,…,加 7=0 m 工ho =6 z=l m 工心=by j J = \,....n Z = 1 s.t. /=! 4布点问题 某市有6个区,每个区都可建消防站,为了肖省开支,市政府希望设置的消防站最少, 但必须保证在该市 任何地区发生火警时,消防车能在15分钟内赶到现场。假定各区的消防站要建的话,就建在区的中心,根据 实地测量, 1区2区3区4区5区6区 1区41016282720 2区10524321710 3区16244122721 4区28321251525 5区27172715314 6区20102125146 请你为该市制左一个设置消防站的最省的讣划。建模并求解。 解:本题实际上是要确泄各个区是否要建立消防站,使其既满足要求,又最卩省。这自然可 引入0?1变量,故设 1,当在第/区建消防站时 (?=]2 6) 0,当不在第/区建消防站时 若1区发生火警,按照“消防车要在15分钟内赶到现场”的要求,贝IJ I 区和2区至少 应有一个消防站,即為+£?1。同理得: x x +x 2 +x 6 > 1, 从而得模型为: 7=1 X, +x 2 >1 x l + x 2 +x 6>\ x 3+x 4>\ S ?f S x 3+x 4+x 5>\ x 4 + x 5 +x 6 > 1 x 2 + x 5 +x 6 > 1 Xj =0,1,( J = 1,2,…,6 丿 若满足第一、三个约束条件,则必然满足第二、四个约朿条件,故后者 从而化简得: tnin f = X X j /-I Xj +x 2 > 1 x 3 +x 4 > 1 s.t.< x 4 + x 5 + x 6 > 1 x 2 + x 5 + x 6 > 1 Xj =0,1,( j = 1,2,…,6 丿 此模型当然可用软件求解,但由于比较简单,故可直接试算。若要求只有一个厂=1,则显 然不可行;若要求只有两个Xj =1 f 则有唯一可行解x 2 =x 4 = 1, Xj = = x 5 = x 6 = 0 , 故这就是最优解,即只需在2区和4区建立消防站。 附程序: c=[1 11111] a=-[1 1 OOOOjOO 1 1 0 O;O 0 0 1 1 1;O 1 00 1 1] b=-[1 1 1 1] v=zeros(1,6) u=[1 11111] X := 6 目标是f = X X J 最少。 以下考虑约束条件。 x 3+x 4> 1, x 3+x 4 +x 5 > 1, x 4+x 5+x 6 > 1, x 2 + x 5 + x 6 >\ 再仔细观察知, 是多余的,可省略。 [x fval]=linprog(c,a,b,[],[],v,u) Optimization terminated? x = 0.0000 1.0000 0.0000 1.0000 0.0000 0.0000 fval = 2.0000 5配套问题 设有门个车间要生产m种产品,第j车间每天生产第i种产品至多Gj件(即全天只安排生产产品,而不安排生产其他产品时的最大产量),假设这m种产品每种一件配成一套, 问如何安排生产任务才能使生产出的成套产品最多? 设Xij车间j安排用于生产产品i的数量,Z——每天生产的成套产品数目, /1 max f = z = min V x. i f>Z s X.. < fl.. x.. > o a = i,…,〃2;/?=i,…,) 模型改进为: max f = z (i = h…") j" Xjj 5 ? A;y > 0 (; = 1,???,〃叮=1,???/) 模型问题:没有将一天的生产时间约束考虑到! 设Xfj ——车间八安排用于生产产品f的时间(占全天的比例) z—一每天生产的成套产品数目 则夠X"——车间/每天生产产品j的数目。例如:车间2每天至多生产某产品6件,若安排1/3天时间去生产,则可产出2件),而工cgxq一一每天全厂产出产品了的总量。y-i 12■ ■ ■ j ■ ■ ■n全厂生产产品1的数址 1■ ■ ■■ ■ ■■ ■ ■ X\n a\n ■ ■ ■ 2 ■ ■ ■■ ■ ■ X2j a2j ■ ■ ■ X2n a2n ■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■ ?1■ ■ ■■ ■ ■■ ■ ■ X A n 为/j全厂毎天生产产品1的总 /?I fit ■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■■ ■ ■m■ ■ ■■ ■ ■■ ■ ■■ ■ ■ 车剛毎天生产产品I的总扯 nt 2>曲车间j i-i 每天生产产品i的 总 max Z ( z = min V ) 土ClijXijNZ = Hl s*龙如兰1(丿=12…力丿/-I x y > O (i = 1,2,???冲,?丿=1,2,.../?丿Z>0 整数 其中常数1表示1天。 注:(H此模型着重考虑安排生产的时间:(2)从实际情况考虑,安排生产的时间必须是每件产品耗用生产时间的整数倍才合适。例如,每3分钟生产一件产品,安排5分钟,也只能生产、件,不能认为生产了\.67件。模型(*)没有考虑到这些因素,故是不合适的。 (2)建模方法(二)一一改进严) 显然,—一一车间/生产每件产品i的耗用时间(天九从以上分析应有 a:: (?是非负整数) 从而令WjfjXij , Wj是非负整数,表示车间J每天生产产品i的数目,将它代入(°)后得max Z 工y nZ仃=12???川丿 冃tJ / 、 " 1 22 —* 1 (j=1,2,…〃丿 /=1偽丿1 y…>0 整数(,=1,2,…川;j = 1,2,.../?丿 " Z>0整数 这是一个整数线性规划模型。 注:此模型着重考虑安排生产产品的数目。 四、数学规划在数学建模中的应用举例 1.钻井布局 勘探部门在某地区找矿。初步勘探时期已零散地在若干位苣上钻井,取得了地质资料。进入系统勘探时期后,要在一个区域内按纵横等距的网格点来布宜井位,进行“撒网式”全而钻探。由于钻一口井的费用很髙,如果新设讣的井位与原有井位重合(或相当接近),便可利用旧井的地质资料,不必打这口新井。因此,应该尽量利用旧井,少打新井,以丹约钻探费用。比如钻一口新井的费用为500万元,利用旧井资料的费用为10万元,则利用一口旧井就节约费用490万元. 设平面上有n个点R,其坐标为(a f,b,),j=1,2,...,n,表示已有的n个井位。新布置的井位是一个正方形网格N的所有结点(所谓“正方形网格”是指每个格子都是正方形的网格:结点是指纵线和横线的交叉点)。假左每个格子的边长(井位的纵横间距)都是1单位(比如100米)。整个网格是可以在平而上任意移动的。若一个已知点Pi与某个网格结点Xi的距离不超过给左误差£(=0.05单位),则认为Pi处的旧井资料可以利用,不必在结点X’ 处打新井。 为进行辅助决策,勘探部门要求我们研究如下问题: 1)假左网格的横向和纵向是固泄的(比如东西向和南北向),并规左两点间的距离为英横向距离(横坐标之 差绝对值)及纵向距离(纵坐标之差绝对值)的最大值。在平而上平行移动网格N,使可利用的旧井数尽可能 大。试提供数值计算方法,并对下而的数值例子用计算机进行计算。 2)在欧氏距离的误差意义下,考虑网格的横向和纵向不固左(可以旋转)的情形,给岀算法及计算结果。 3)如果有门口旧井,给出判左这些井均可利用的条件和算法(你可以任意选左一种距离)。 数值例子n“2个点的坐标如下表所示: I123456789101112 ②问题分析 布局问题是在一定约束条件下的最优化问题,勘探部门要求尽可能多地利用旧井,因此我们必须围绕这个问题来移动网格。 问题(1):网格的横向和纵向是固定的,网格只能平移。如果考虑网格上所有的结点通过平移后与旧井位的距离,由于网格结点数较多,运算量较大,虽然可以解决问题,但并非一个好的解法。因此我们从旧井位考虑,通过四舍五入法找到与其相邻的网格结点(相互关系如图所示),然后我们通过两个控制变量(角度与距离)对这些结点进行平移,使网格结点尽可能多的与旧井位重合。 问题(2):网格不固左移动,需要我们考虑平移和旋转,对于旋转我们同样从旧井位出发,通过旧井位与网格结点的相对关系,不难得出:旧井位的旋转是网格旋转的逆过程, 因此我们利用旧井位旋转来替代网格的旋转问题,然后,根据问题(1)的方法,再解决网格平移问题,这样就将网格的移动问题简化成为旧井位的旋转问题和网格的平移问题。 问题(3):以上两个问题都是在正方形网格情况下进行的求解,在此问中我们不妨假设网格是长方形的,通过调整长和宽的比例关系使n个旧井位都可以被利用。因此问题(3), 就是在问题(2)的基础上,使n 个旧井位都可以被利用的长方形的长、宽比例问题。具体实现方法如下:把目标函数值定为12,通过对问题(2)的反向求解来确定长方形的长、宽比例。 ③模型假设和符号说明 模型假设: 1. 若一个旧井位与某个网格结点不超过0.05单位,则视为重合: 2. 网格是足够大且平铺的: 3. 模型采用直角坐标系,在网格未移动时,网格的横向和纵向就是直角坐标系的两坐标轴方向。 4. 逆时针旋转为正,顺时针方向为负: 5. 长方形网格的长是4单位。 符号说明: i...i 旧井1=1,2, (12) Pi...第i个旧井位; ai,bi...第i个旧井位的横坐标、纵坐标; Ai,Bi...第i个旧井位旋转的横坐标、纵坐标; Wi...与第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)决泄,即而由要求知 1-1 fl 0< d(i) < 0.05 ti(d)=|o沏皿⑴ d(i)的表达式在下面给出。 由于问题<1)是问题(2)的特殊情况,所以这里只给岀问题(2)的模型。我们先将旧井位进行旋转。 第i个旧井位到坐标原点的距离为 p(i)= Jo" +的)2 . (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⑴=y(i)+s*sin(an) (6) 所以,网格结点与旧井位的渓差距离为 d(i)= ⑺ 综上所述,本问题的数学模型可表示为: 目标函数:maxT=工〃(〃) /-] 约束条件是满足(/)-(7)式。 注:以上模型说明:当bn=0, d(i)=max{|xl(i)-A(i)|, |yl(i)-B(i)|}时,目标函数的解即为问题(1)的解。 ⑤模型的求解 我们的目标是使目标函数最大化,这是一个有约朿的优化问题,我们采用穷举法来求它的极值。所谓穷举法就是逐点确定寻优方向,然后利用固左的步长,进行搜索的方法。为使目标函数值最大,我们列出主要步骤如下: 1. 泄一个能包含所有解的初始范用与固左的搜索步; 2. 据运行时间,再对搜索步长用穷举法进行精度调整。 根据这种方法我们得到了问题答案: 问题(1):最多可利用的旧井位数t=4,网格的平移方向an0=5.400000(弧度),网格的平移距离 sO=O.583OOO?位); 问题(2):最多可利用的旧井位数U6,网格的平移方向an0=2.760000(弧度),网格的平移距离 sO=O.23OOOO(单位),网格旋转角度:bnO=O.78OOOO(弧度)(旧井点的相对旋转角度为: 5.500000弧度)。 ⑥模型的改进 降落伞的选择 一.问题的重述: 为向灾区空投一批救灾物资共2000kg,需选购一些降落伞。已知空投高度为500m,要求降落伞落地时的速度不能超过20m/s,降落伞的伞面是半径为r的半球而,用每根长L 共16根绳索连接的重m位于球心正下方球而处,如下图: 每个降落伞的价格由三部分组成,伞而费用G由伞的半径厂决泄,见下表:绳索费用C2由绳索总长度及单价4元每。 降落伞在降落过程中除受到重力外,还受到空气的阻力,可以认为与降落的速度和伞的而积的乘积成正比。为了确 试确立降落伞的选购方案,即共需多少个伞,每个伞的半径多大(在给左半径的伞中选), 在满足空投要求的条件下,使费用最低。 二、问题的分析: 根据题意,每种伞的价格是确定的。要确左降落伞的选购方案,即共需多少个伞,每个伞的半径多大(在给立半径的伞中选),在满足空投要求的条件下,使费用最低。首先,必须知道每种伞在满足空投要求的条件下最大的载重SAI(r),然后就是一个线性整数规划了。欲得到Al(r).必须先求出空气阻力系数k,然后根据运动方程得出A1(r)o最后运用分支泄界法求解线性整数规划,得出问题要求的解匚 三、基本假设: 1.救灾物资2000kg可任意分割: 2?降落伞落地时的速度不能超过20m/s; 3. 降落伞与绳索的质量可以忽略: 4. 降落伞下落过程中,只受到重力和空气阻力的作用: 5. 空气阻力的阻力系数k是世值,与苴他因素无关。 四、符号说明: A1(r)…… 半径为r的伞在满足空投要求的条件下最大的载重量: t……降落伞从开始下降开始记时的时间; k……空气阻力系数 H(r)……降落伞从降落位置到t时刻所下降的距离; m.... 降落伞负重重呈:: g……重力加速度; s-降落伞伞而面积; n r……选购的半径为r的降落伞的个数 五、模型的建立与求解 计算每种伞的单价如下:单位为元,C2=^2rxl6x4 2.定空气阻力系数 22 作岀的关系图 SCO 450 4C0 393 3C0 2£0 2C0 ISO ICO EO 图1 从图1可看岀:一方而,H(m)~t(s)在后阶段基本是线性关系,即降落伞作匀速运动。 ms=kVs 有又由于若在500米的物体做匀速运动,则v = - = —= Um/s , t 30 代回mg=kVs中,估算岀后29 另一方面,根据题意,物体在降落过程中一直做匀速直线运动是不可能的。 题中注意到:降落伞在下降过程中只受到重力和空气阻力的作用,而且初速度为0, 由牛顿第二泄理知:l dV(t) . kvs F = ma = m ------- = IHQ_ kvs ■ => --------- = g --------- dt dt m (dV(t) kvs ⑴ V(0) = 0 kst (2) 解得:%)=竺一曲i ks ks 由(2)式,并代入k=2.9, r=3m, m=300kg,取g=9.8m/s2> s=2 Jt r2,作出下图2 图2 从图2可以看出当0 快就接近极限值mg/ks.当宀9时,一竺二= 竺,所以9秒以后 ks;甞ks kse m 可以看作近似的匀速运动。下而就9秒以后的数据运用最小二乘法进行线性拟合。 设H(t)=Vt+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)> 并由ms=kVs得到k=2?9575。 3、瞬时速度与高度 设降落伞从降落位置到t时刻所下降的距离为H (t),则有 4.求半径为r的降落伞在满足空投要求的条件下最大的载重量M(r). 叩)与/直接相关,也与加相关。加越大,则v 越大。由(2)式知:V(m)是关于m 的增函 数。当v 最大时,川也达到最大M ,即为最大的载重量。 特别的,在给左从500m 的高空空投时,降落伞落地瞬间的速度在给左g, k, s 后又有等式 约束H(t)=500,即 z 〃妙,加ge ,n 加g ks k= k = 的情况下,V(m)是关于m 的增函数:反之,其反函数也是关于V 的增函数。所以要求半径 为r 的降落伞在满足空投要求的条件下最大的载重疑Al(r),就是要在V 取最大值时取得, 即取V=20m/s 时求岀指定半径r 的M(r),于是得到如下方程组: ks ks kxi (5) 由此导出: ”、 厶 ksV x /zf -> 八 mV ,八 H(f)= 一〃厂g ? ln(l ----- )/伙飞~)一 ----- (6) mg ks 如前所述,取H(t)=500, V=20m/s,得到方程: L s y 500 = 一加2g ? ln(l- —)/(k 2s 2)- mg 将参数g 二9?& k=2.9575代入上式,有 由 s=2 n r 2,分别代入 r=2t r=2.5: r=3: r=3.5; r=4,并调用 Matlab 的命令 solve> 分别解 得半径为r M(2)kg M(2?5)kg M(3)kg M(3?5)kg M(4)kg 151.6947 237.0229 341.3130 464.5649 606.7787 (取整)151 237 341 464 606 3- 4.原问题就是如下的一个线性整数规划问题: 设nr 为选购的半径为r 的降落伞的个数,则有 ksi (7) 500 =-9.8-nz 2 ln(l- 20?2.9575s 9.8/n )/(2.95752 ? 一 20m 2.9575-5 ks t nin{ 446//2 +596/?25 +822n3+1177 n35 +1562 — } 151弘+ 237小s + 341? + 464/:3 s+ 606/“ > 2000 s.t. < ?「(8) n r,r = 2,2.5,3,3.5,4 n r e N 运用分支定界法可求得:n2 = n2S= “35 = n4 = 0, n3 = 6 总费用为4932元,即在需要空投2000kg的情况下,需要选购半径为3米的降落伞6把。 六.模型的检验及推广 1.在求解空气阻力系数k时,我们分析数据得出在运动后期降落伞作近似的匀速运动,并由此为前提对数据进行 可以看出,M/s几乎是常数,又有m S=kVs是降落伞后期运动为匀速运动的充分必要条件,即m/s=kV/g为常数,而k, g为常数,所以降落伞在运动后期为近似匀速运动。因此,在求解空气阻力系数k时,假设后期运动为近似匀速运动是合理的。 2. 从图2可以看岀降落伞的大幅度加速过程很快就结束了,对给定的伞和给泄的承载质量很快就进入近似匀速运动,而且速度与空投高度基本上是无关的,所以空投高度并不十分重要,只要能保证空投位苣准确,高一些投放也可以,这样可降低空投髙度。 3. 题目中要求的是空投2000kg的救灾物资,若数据有变动,例如3000kg, 5000kg?则只需将(8)式的第一个约束右端改为3000kg, 5000kg,然后再运用分支左界法可求得结果。