一、线性规划问题 题目1：
评分：
一车队有8辆车，这8辆车存放在不同地点，队长要派其中5辆到个不同的工地去运货。
各车从存放处调到装货地点所需费用列于表2.5.问应选用哪五辆车调到何处去运货，才能
使各车从所在地点调到装货地点所需的费用最少（要求分别用matlab和Lingo编程求解）？
表2.5 调车费用表-
（1）
解：记cij表示第j号车调到装货地点i所需的费用。引进0-1变量，
?1,第j号车调到装货地点ixij??
0，第j号车厢没有调到装货地点?
建立如下的0-1整数规划模型：
5
8
min??cijxij
i?1j?1
?5
??xij?1,j?1...8?i?1?8
s.t??xij?1,i?1...5 ?j?1
?x?0或1,i?1,2...5;j?1,2....8?ij?
（1）编写的matlab程序如下 clc,clear
c=[30 25 18 32 27 19 22 26
29 31 19 18 21 20 30 19 28 29 30 19 19 22 23 26 29 30 19 24 25 19 18 21 21 20 18 17 16 14 16 18]; c=c(:);
a=zeros(8,40);
for j=1:8
a(j,[(j-1)*5+1:5*j])=1; end
b=ones(8,1);
d=zeros(5,40);
for i=1:5
d(i,[i:5:40])=1;
end
e=ones(5,1);
[x,fval]=bintprog(c,a,b,d,e); x=reshape(x,[5,8]) fval 运行结果如下：
x =
(1,3) 1
(2,4) 1
(3,5) 1
(5,6) 1
(4,7) 1
fval =
87
所以最优方案为x13?x24?x35?x56?x47?1,最优值为 87。
Lingo程序如下
model:
sets:
col/1..5/:i;
var/1..8/:j;
links(col,var):c,x;
endsets
data:
c=30 25 18 32 27 19 22 26 29 31 19 18 21 20 30 19 28 29 30 19 19 22 23 26
29 30 19 24 25 19 18 21
21 20 18 17 16 14 16 18;
enddata
min=@sum(links:c*x);
@for(col(i):@sum(var(j):x(i,j))=1);
@for(var(j):@sum(col(i):x(i,j))<=1);
@for(links:@bin(x));
end
运行结果为
Objective value: 87.00000
X( 1, 1) 0.000000 30.00000 X( 1, 2) 0.000000 25.00000 X( 1, 3) 1.000000 18.00000 X( 1, 4) 0.000000 32.00000 X( 1, 5) 0.000000 27.00000 X( 1, 6) 0.000000 19.00000 X( 1, 7) 0.000000 22.00000 X( 1, 8) 0.000000 26.00000 X( 2, 1) 0.000000 29.00000 X( 2, 2) 0.000000 31.00000 X( 2, 3) 0.000000 19.00000 X( 2, 4) 1.000000 18.00000 X( 2, 5) 0.000000 21.00000 X( 2, 6) 0.000000 20.00000 X( 2, 7) 0.000000 30.00000 X( 2, 8) 0.000000 19.00000 X( 3, 1) 0.000000 28.00000 X( 3, 2) 0.000000 29.00000 X( 3, 3) 0.000000 30.00000 X( 3, 4) 0.000000 19.00000 X( 3, 5) 1.000000 19.00000 X( 3, 6) 0.000000 22.00000 X( 3, 7) 0.000000 23.00000 X( 3, 8) 0.000000 26.00000 X( 4, 1) 0.000000 29.00000 X( 4, 2) 0.000000 30.00000 X( 4, 3) 0.000000 19.00000 X( 4, 4) 0.000000 24.00000 X( 4, 5) 0.000000 25.00000 X( 4, 6) 0.000000 19.00000 X( 4, 7) 1.000000 18.00000 X( 4, 8) 0.000000 21.00000 X( 5, 1) 0.000000 21.00000
X( 5, 2) 0.000000 20.00000 X( 5, 3) 0.000000 18.00000 X( 5, 4) 0.000000 17.00000 X( 5, 5) 0.000000 16.00000 X( 5, 6) 1.000000 14.00000 X( 5, 7) 0.000000 16.00000 X( 5, 8) 0.000000 18.00000
所以最优方案为x13?x24?x35?x56?x47?1,最优值为 87。
第2题
某单位需要加工制作100套钢架，每套用长为2.9m，2.1m和1m的圆钢各一根。已知原料长6.9m,问应如何下料，使用的原材料最省。 解：列出所有可行性方案如下表
第一种方案x1：切割2.9米一截,2.1米一截,1米一截，浪费0.9米； 第二种方案x2:切割2.9米一截，1米4截，浪费0米； 第三中方案x3：切割2.9米2截,1米1截，浪费0.1米； 第四种方案x4：切割2.1米3截，浪费0.6米；
第五种方案x5:切割2.1米2截，1米2截，浪费0.7米； 第六种方案x6：切割2.1米1截，1米4截，浪费0.8米； 第七种方案x7：切割1米6截，浪费0.9; 设：使用第i方案xi次（i=,3,4,5,6,7） Min=x1+x2+x3+x4+x5+x6+x7；
?x1?x2?2x3?100
?x1?3x4?2x5?x6?100?
s.t?
?x1?4x2?x3?2x5?4x6?6x7?100??x1,x2,x3,x4,x5,x6,x7?0
model: sets:
row/1..3/:b; col/1..7/:c,x; links(row,col):a;
endsets
data:
c=1 1 1 1 1 1 1;
a=-1 -1 -2 0 0 0 0
-1 0 0 -3 -2 -1 0
-1 -4 -1 0 -2 -4 -6;
b=-100 -100 -100;
enddata
min=@sum(col:c*x);
@for(row(i):@sum(col(j):a(i,j)*x(j))<=b(i));
@for(col(j):@gin(x));
end
运行结果：
Objective value: 91.00000
X( 1) 0.000000 1.000000 X( 2) 12.00000 1.000000 X( 3) 44.00000 1.000000 X( 4) 33.00000 1.000000 X( 5) 1.000000 1.000000 X( 6) 0.000000 1.000000 X( 7) 1.000000 1.000000
则最少根数为91根，选择的方案为第二种方案12根，第三种方案44根，第四种方案33根，第五种方案1根，第七种方案1根。
二、图与网络问题
评分：
第一题
已知有6个村庄，各村庄的小学生人数如表4.11所列，各村庄间的距离如表4.5所示。现在计划建造一所医院和一所小学，问医院应建在哪个村庄才能使最远村庄的人到医院看病所走的路最短？又问小学建在哪个村庄使得所有小学生走的总路程最短？
解：（1）建立赋权完全图G=（V，E，W），其中V=（），权重邻接矩阵为：
?027?????20468?????74013?? w???
??61016???83103??????630????
利用Floyd算法可以求出任意两个村庄的最短距离矩阵：
11?? 0 2 6 7 8 ? 2 0 4 5 6 9???? 6 4 0 1 2 5? A?(aij)6?6?? 7 5 1 0 1 4??? 8 6 2 1 0 3???11 9 5 4 3 0????
其中表示第i村庄到第j村庄的最短距离。
A的第j列表示其他各村庄到该村的最短距离，第j列最大值最远的村到该村的最短距离，由矩阵su可知，每一列的最大值为11,9,6,7,8,11；所以医院应建在第三村庄使得最远村庄的人到医院看病所走的路最短，为6. 相应的matlab程序如下： clear;clc;
a=zeros(6);a(1,2)=2;a(1,3)=7;a(2,3)=4;a(2,4)=6;a(2,5)=8;a(3,4)=1;a(3,5)=3;a(4,5)=1;a(4,6)=6;a(5,6)=3; a=a';a=sparse(a);
A=graphallshortestpaths(a,'Directed',0) su=max(A)
（2）A矩阵的第i行表示其他村庄到第i村庄的最短距离，设表示第i村庄的人数（i=1,2..6）,若在第i村上学，则所有小学生走的总路程为，则当时，即为所求。
用matlab求解得：
s=[ 0 2 6 7 8 11; 2 0 4 5 6 9; 6 4 0 1 2 5; 7 5 1 0 1 4; 8 6 2 1 0 3; 11 9 5 4 3 0]; a=[50 40 60 20 70 90]; s=A*a'
结果如下列表：
第二题
乙两个煤矿分别生产煤500万吨，供应A、B、
C三个电厂发电需要，各电厂用量分别 为300、300、400（单位:万吨）。已知煤矿之间、煤矿与电厂之间以及各电厂之间相互距离（单位：km）如表4.13～表4.15所列。煤可以直接运达，也可以间接运达试确定从煤矿到各电厂间煤的最优调配方案。
如图所示把甲，乙，A，B，C分别用1，2，3，4，5表示。由题意知这是一个多源多汇图。我们把源点聚为s点，收点聚为t点。这题就形成了求从s点到t点的最大流最小费用问题。 可得到如下数学模型：
model:
sets:
nodes/s,1,2,3,4,5,t/:d;
arcs(nodes,nodes)/s 1,s 2,1 2,1 3,1 4,1 5,2 1,2 3,2 4,2 5,3 4,3 5,4 3,4 5,5 3,5 4,3 t,4 t,5 t/:b,c,f;
endsets
data:
d=1000 0 0 0 0 0 -1000;
b=0 0 120 150 120 80 100 60 160 40 70 100 50 120 100 150 0 0 0; c=500 500 10000 10000 10000 10000 10000 10000 10000 10000 10000 10000 10000 10000 10000 10000 300 300 400;
enddata
min=@sum(arcs:b*f);
@for(nodes(i):@sum(arcs(i,j):f(i,j))-@sum(arcs(j,i):f(j,i))=d(i));
@for(arcs:@bnd(0,f,c));
end
求解结果：
F（1,4）=300，F（1,5）=200，F（2,3）=300，F(2,5)=200，甲给B，C各运300,200万吨，乙给A,C各运300，200万吨。1000万吨媒运输的总路程为78000（Km*万吨）。
三、微分方程问题
第一题
根据经验，当一种新商品投入市场后，随着人们对它的拥有量的增加，其销售量s(t)下降的速度与s(t)成正比。广告宣传可给销量添加一个增长速度，它与广告费a(t)成正比，但广告只能影响这种商品在市场上尚未饱和的部分（设饱和量为M）。建立一个销售s(t)的模型。若广告宣传只进行有限时间T，且广告费为常数a，问s(t)如何变化？
解：设?为s(t)的衰退因子，p为a(t)对s(t)的影响系数。又因为s(t)的下降速度与s（t）是成正比的。所以在没有广告的影响下评分： ds???s(t)，又因为广告费可以给销量增添一dt个增长速度，与a（t）成正比，且只影响未饱和部分。
所以建立如下模型：
dss(t)?pa(t)(1?)??s(t)，（1） dtM
（1）当s(t)=M或a(t)=0时：
ds???s(t) dt
（2）假设采用如下的广告策略：
?a/r,0?t?r a(t)??0,t?r?
当0?t?r时，将a(t)=a/r带入（1）式得：
dspas(t)?(1?)??s(t)可推出 dtrMdspapa?(??)s(t)? dtrMr
papa??=b, 令=c. rrM
原式可化为如下的常微分方程：
ds?bs(t)?c dt
若令s（t）=s0，求出解为：用matlab求
y=dsolve('Ds+b*s=c','s(0)=s0','t')解得：
y =(c - (c - b*s0)/exp(b*t))/b； 简化为：s(t)?c(1?e?bt)?b0s0e?bt b
解为：
c?0?t?r ?s(t)?(1?e?bt)?b0s0e?bt，?b?(r?t)?s(t)?s(r)e,t?r?
第二题
2有高为4m的半椭球形容器，水从它的底部小孔流出。小孔横截面积为1cm。开始时容器内
盛满了水，求水从小孔流出过程中容器里水面的高度h（水面与孔口中心的距离）随时间 t
变化的规律。
解：以底部中心为坐标原点，垂直向上为坐标的正向建立坐标系。由能量守恒原理得到：
mgh?12mv, 2
解得v?2gh。
设在微笑时间间隔[t,t?dt]内，水面高度由高度h降到h+dh（这里dh为负值）， 由物质守恒原理得：
是s(?dh)?s0vdt （1），
其中：s为水面所在的截面面积，s0为底部小孔的面积。有：
(4?h)2
s?2?(1?)（2） 16
把v?2gh，s0?0.0001和（2）代入（1）式，化简得：
(4?h)2
20000(1?)dt?(?dh) 2gh
在考虑初始条件，得到如下微分方程模型：
?(4?h)2
20000?(?1)??dt ??dh2gh??t(1)?0?
解得微分方程的解：
t=dsolve('Dt=20000*pi*((4-h)^2/16-1)/sqrt(2*g*h)','t(4)=0','h'); t=simple(t),pretty(t)
53解得：t?h2?40h2?244)
四、目标规划问题 评分：
一．最近，某节能灯具厂接到了订购16000套A型和B型节能灯具的订购合同，合同中没有对这两种灯具各自的数量做要求，但合同要求工厂在一周内完成生产任务并交货。根据该厂的生产能力，一周内可以利用的生产时间为20000min，可利用的包装时间为36000min。生产完成和包装完成一套A型节能灯具各需要2min；生产完成和包装完成一套B型节能灯具分别需要1min和3min。每套A型节能灯具成本为7元，销售价格为15元，即利润为8元；每套B型节能灯具成本为14元，售价为20元，即利润为6元。厂长首先要求必须要按照合同完成订货任务，并且既不要有不足量，也不要有超过量。其次要求满意的销售额尽量达到或接近275000元。最后要求在生产总时间和包装总时间上可以增加，但超过量尽量地小。同时注意到增加生产时间要比增加包装时间困难得多。试为该节能灯具厂制定生产计划。
解：根据题意，确定目标和优先级：
第一级目标：首先要求必须要按照合同完成订货任务，赋予优先因子为p1;
第二级目标：其次要求满意的销售额尽量达到或接近275000元，赋予优先因子为p2； 第三级目标：最后要求在生产总时间和包装总时间上可以增加，但超过量尽量地小，赋予优先因子为p3；
建立相关目标约束，设生产A型和B型节能灯各x1，x2套。
（1）订单数量约束。
设d1和d1分别表示未达到和超过订单量的偏差量，又因为只要求完成16000套交货数量即可，则数量目标约束为 ??
????min{d1?d1} ?????x1?x2?d1?d1?16000
（2）销售利润目标。
因为要求利润尽量完成和接近275000元，因此得
???min{d2} ?????15x1?20x2?d2?d2?275000
（3）生产总时间和包装总时间的限制。
因为要求在生产总时间和包装总时间上可以增加，但超过量尽量地小，同时注意到增加生产时间要比增加包装时间困难得多，则可设二者的加权系数为0.4和0.6，则
??min{0.4d3??0.6d4}????2x1?x2?d3?d3?20000
???2x?3x?d?d?360001244?
写出目标规划的数学模型，即
???minz?p1(d1??d1?)?p2d2?p3(0.4d3?0.6d4)
?x1?x2?d1??d1??16000????15x1?20x2?d2?d2?27500? s.t?2x1?x2?d3??d3??20000
???2x?3x?d?d?360001244??x,x,d?,d??0,i?1,2,3,4.?12ii
写出相应的lingo程序如下：
sets:
a/1..3/:p,z,goal;
b/1..2/:x;
c/1..4/:g,dplus,dminus;
d(c,b):e;
obj(a,c)/1 1,2 2,3 3,3 4/:wplus,wminus;
endsets
data:
ctr=?;
goal=？？0;
g=16000 275000 20000 36000;
e=1 1 15 20 2 1 2 3 ;
wplus=1 0 0.4 0.6;
wminus=1 1 0 0;
enddata
min=@sum(a:p*z);
p(ctr)=1;
@for(a(i)|i#ne#ctr:p(i)=0);
@for(a(i):z(i)=@sum(obj(i,j):wplus(i,j)*dplus(j)+wminus(i,j)*dminus(j)));
@for(c(i):@sum(b(j):e(i,j)*x(j))+dminus(i)-dplus(i)=g(i)); @for(a(i)|i#lt#@size(a):@bnd(0,z(i),goal)); end
运行结果如下：
Objective value:3800.000
X(1) 9000.000 0.000000 X(2) 7000.000 0.000000 DPLUS(3) 5000.000 0.000000 DPLUS(4) 3000.000 0.000000
求得的最优解为生产A型灯具9000套，B型灯具7000套。生产时间需要增加5000min，包装时间需要增加3000min，即可完成任务。
五、评价与预测问题
一．某商品的生产需要甲乙两种原料，产品利润，以及甲乙两种原料的市场供给等数据如表15-7所列，试预测20xx年，甲的供应量为400KG，乙的供应量为500KG时的产品利润（要求建立灰色GM（1，1））。
表15-7 原始数据表
评分：
量，分别表示利润，甲原料，乙原料。对应三个原始数据，设参数数据分别为
00000x1?(x1(1),x1(2),x1(3),x1(4),x10(5)) 000000x2?(x2(1),x2(2),x2(3),x2(4),x2(5))
000000x3?(x3(1),x3(2),x3(3),x3(4),x3(5))
（1） 求级比?i(k)（i=1,2,3），有
x(0)(k?1) ?(k)?，
x(0)(k)
可以得?1?(?1(2),?1(3),?1(4),?1(5))=（0.5748,0.7262,0.9279,0.6351）
?2?(?2(2),?2(3),?2(4),?2(5))=（0.6336,0.7279,0.9231,0.6373）
（1）级比判断。由于?i(k)?（0.7165,1.3307），（i=1,2,3），所以需要对序列xi0做必要的变换处理，使其落入可覆盖内。即取适当的常数c，作平移变换
00对x1作平移变换，y1(k)?x10(k)?6000，k=1,2….5，
'则?1=（0.7621，0.8258，0.9529，0.7270）?（0.7165,1.3307）
000对x2作平移变换，y2(k)?x2(k)+100，k=1,2….5,
则?'
2=（0.7922，0.8250，0.9492，0.7266）?（0.7165,1.3307）
000对x3作平移变换，y3(k)?x3(k)+150，k=1,2….5,
?'
3=（0.8177，0.9452，0.9364，0.7383）?（0.7165,1.3307）
由于所有的?;
i(k)?（0.7165,1.3307），k=2,3…5,故可以用xi0（i=1,2,3）作令人满意的GM（1,1）
建模。
2）GM(1,1)建模
（1）对平移变换后数据xi0（i=1,2,3）作一次累加，得到
x1?1?=（10383，13625，16500，17316，23818）
?1?=（183,414,694,989,1395） x2
?1?=（296,658,1041,1450,2004） x3
（2） 构造数据矩阵Bi（i=1,2,3）及数据向量Yi（i=1,2,3），有 因为x1的均值生成序列为z1
（12004,15062.5,16908,20567）
?1??1???z1?1?(2)x2(2)x3(2)???12004414??1????1??1?.5694?z(3)x(3)x(3)123?=??15062所以B???1??1??1???z1(4)x2(4)x3(4)???16908989??1????1??1??z(5)x(5)x(5)??23?1???205671395?1??1?=（z1(2)，z1(3)，z1(4)，z1(5))=?1??1??1??1?658?1041?? 1450??2004?
?y1?0?(2)??13625???0????16500y(3)1? Y???0??=??y1(4)?? 17316???0????y(5) 23818???1???
? -1.9868??a?
?1T????T所以u?b1=?BB?BY? 118.8182 ????????-91.5799??b2??
于是得a= -1.9868，b1=118.8182，b2= -91.5799。
建立模型：
dx(1)?1??1??1.9868x1?1??118.8182x2?91.5799x3 dt
，
及近似时间响应式
x1?1?(k?1)?(x1?0?(1)?
?1?b1?1?b?1?b?1?b?1?x2(k?1)?2x3(k?1))e?ak?1x2(k?1)?2x3(k?1)=(aaaa1.9868k?1?4383+59.8x2(k?1)-46.1x3(k?1))e?1?-59.8x2(k?1)+46.1x3(k?1) ?1?
clc,clear
format long g
x0=[10383 13625 16500 17316 23818
183 231 280 295 406
296 362 383 409 554];
[m,n]=size(x0);
x1_d=cumsum(x0,2)
x11=x1_d(1,:)
z11=0.5*(x11(1:end-1)+x11(2:end))
b=[-z11',x1_d(2,2:end)',x1_d(3,2:end)']
y=x0(1,2:end)'
u=b\y
x1=dsolve('Dx1+a*x1=b1*x2+b2*x3','x1(0)=x10')
x1=subs(x1,{'a','b1','b2','x10'},{u(1),u(2),u(3),x0(1,1)}) x1_s=vpa(x1,7),x1_s=simple(x1_s)
x20=[x0(2,:),400];
x30=[x0(3,:),500];
x21=cumsum(x20);
x31=cumsum(x30);
x1=subs(x1,{'t','x2','x3'},{[0:n],x21,x31})
x10hat=[x1(1),diff(x1)]
epsilon=(x0(1,:)-x10hat(1:end-1))./x0(1,:) delta=abs(epsilon./x0(1,:))
xhat=x10hat(end)
由运行结果可知：预测值为23664.6。
残差r=（0，0.13，-0.08，-0.03，-0.008）
相对误差r’=(0, 9.56e-6,4.89e-6,1.7e-6, 3.39e-7)