数学建模数学规划模型.ppt

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,数学建模数学规划模型,第一节 线性规划模型,（1）目标函数是决策变量的线性函数。,（2）约束条件都是决策变量的线性等式或不等式。,MATLAB,命令命令输入格式及线性规划模型如下：,其中：,x0,是算法迭代的初始点；,nEq,表示等式约束的个,数。,三、,建模举例,营养配餐问题,每种蔬菜含有的营养素成份是不同的，从医学上知道，每人每周对每种营养成分的最低需求量。某医院营养室在制定下一周菜单时，需要确定表,6-1,中所列六种蔬菜的供应量，以便使费用最小而又能满足营养素等其它方面的要求。规定白菜的供应一周内不多于,20,千克，其它蔬菜的供应在一周内不多于,40,千克，每周共需供应,140,千克蔬菜，为了使费用最小又满足营养素等其它方面的要求，问在下一周内应当供应每种蔬菜各多少千克？,表,2-3,问题分析与模型建立,设 分别表示下一周内应当供应的青,豆、胡萝卜、菜花、白菜、甜菜及土豆的量，费用目标函数为：,约束条件：,铁的需求量至少,6,个单位数：,磷的需求量至少,25,个单位数：,维生素A的需求量至少17500个单位：,维生素C的需求量至少245个单位：,烟酸的需求量至少5个单位数：,每周需供应140千克蔬菜，即,0,x,1,40 0,x,2,40 0,x,3,40,0,x,4,20 0,x,5,40 0,x,6,40,问题是满足营养素要求的条件下，所需费用最小，是一个线性规划模型。,利用Matlab软件编程序：,%营养配餐ch21%文件名：ch21 mc=5;5;8;2;6;3;A=(-1)*1,1,1,1,1,1;,0.45,0.45,1.05,0.40,0.50,0.50;,10,28,59,25,22,75;,415,9065,2550,75,15,235;,8,3,53,27,5,8;,0.30,0.35,0.60,0.15,0.25,0.80;,b=(-1)*140;6;25;17500;245;5;,xLB=zeros(6,1);,xUB=40;40;40;20;40;40;,nEq=1;,x0=0*ones(6,1);,x=lp(c,A,b,xLB,xUB,x0,nEq);,disp(,青豆需要的份数,),x(1),disp(,胡罗卜需要的份数,),x(2),disp(,菜花需要的份数,),x(3),disp(,白菜需要的份数,),x(4),disp(,甜菜需要的份数,),x(5),disp(土豆需要的份数),x(6),执行后输出,青豆需要的份数,ans=40,胡罗卜需要的份数,ans=40.0000,菜花需要的份数,ans=,0,白菜需要的份数,ans=,20.0000,甜菜需要的份数,ans=,0,土豆需要的份数,ans=,40,最小费用,ans=,560.0000,背景：,0-1,规划是数学规划的组成部分，起始,20,世纪,30,年代末，七八十年代是数学规划飞速发展时期，无论是从理论上还是算法方面都得到了进一步完善。时至今日数学规划仍然是运筹学领域中热点研究问题，从国内外的数学建模竞赛的试题中看，有近,1/2,的问题可用数学规划进行求解。其中利用,0-1,规划及,0-1,型变量的数学建模问题也为数不少，如,98,年的,投资的收益和风险,，,2004,年的,DVD,在线租赁,等问题，下面我们就来学习,0-1,规划,0-1,型变量在数学建模中的应用。,2.2 0-1规划,0-1型变量,在数学建模中的应用,1、0-1规划,数学规划模型的一般表达式：,整数规划中决策变量只取0或1的特殊情况是0-1规,划。下面通过几个例子说明0-1规划在实际问题中的,应用,。,例2.1 背包问题,有几件物品，编号为 1，2，n。,第 件重为 kg，价值为 元。今有一位装包者欲,将这些物品装入一包，其质量不能超过 kg，问应,装入哪几件价值最大？,解,引入变量,，将 物品装包,，不将 物品装包,于是得问题的模型为,取0或1，i=1，2，n,背包问题看似简单，但应用很广，例如某些投资问,题即可归入背包问题模型。此类问题可以描述为：,投资问题：,设有总额为 元的资金，投资几项事,业，第 项事业需投资 元，利润为 元，问应选,择哪些项投资总利润为最大？,例2.2,某钻井队要从以下10个可供选择的井位中确,定5个钻井探油，使总的钻探费用为最小。若10个,井位的代号为 ，相应的钻探费用为,，并且井位选择要满足下列限制条件：,（1）或选 和 ，或选 ；,（2）选择了 或 就不能选 ，反之亦然；,（3）在 中最多只能选两个。,试建立其数学模型：,解,引入变量,选择,不选择,于是以上问题的数学模型为：,投资的收益和风险,这是,1998,年全国大学生数学建模竞赛的,A,题,问题如下：市场上有,n,种资产（股票、债券、,）,S,i,（,i,=1,，,，,n,）供投资者选择，某公司有数额为,M,的一笔相当大的资金可用作一个时间的投资。,公司财务分析人员对这,n,种资产进行了评估，估,算出在这一时期内购买,S,i,有平均收益率为,r,i,，并预,测出购买,S,i,的风险损失率为,q,i,。考虑到投资越分散,总的风险越小，公司确定，当用这笔资金购买若,干种资产时，总体风险可用所投资的,S,i,中最大的一,个风险来度量。购买,S,i,要付交易费，费率为,p,i,，并且,当购买额不超过给定值,u,i,时，交易费按购买,u,i,计算,（不买当然无须付费）。另外，假定同期银行存款,利率是,r,0,，且既无交易费又无风险（,r,0,=5%,）。（,1,）已知,n,=4,时的相关数据如下：,试给该公司设计一种投资组合方案即用给定的资金,M,，有选择地购买若干种资产或存银行生息，使净,收益尽可能大，而总体风险尽可能小。,（,2,）试就一般情况对以上问题进行讨论，利用,以下数据进行计算。,在,AD,边等距地设置,7,个波源,R,i,（,i,=1,，,，,7,），,BC,边对等地安放,7,个接收器,S,j,（,j,=1,，,，,7,），记录由,R,i,发出的弹性波到达,S,j,的时间,ij,秒）,，,见表,2-3,。,2,、,0-1,型变量在数学建模中的应用,1,、空洞探测问题,1.1,问题的提出,这是,2000,年全国大学生数学建模竞赛的,D,题。山体、隧洞、坝体等的某些内部结构可用弹性波测量来确定。一个简化问题可描述为，一块均匀介质构成的矩形平板内有一些充满空气的空洞，在平板的两个邻边分别等距地设置若干波源，在它们的对边对等地安放同样多的接收器，记录弹性波由每个波源到达对边上每个接收器的时间，根据弹性波在介质中和空气中不同的传播速度，来确定板内空洞的位置。现考察如下的具体问题：,一块,240,（米）,240,（米）的平板,（如图,21,）,在AB边等距地设置7个波源P,i,（,i,=1，7），CD边对等地安放7个接收器Q,j,（,j,=1，7），记录由P,i,发出的弹性波到达Q,j,的时间,t,ij,（秒）,，,见表,2-2,；,在,AD,边等距地设置,7,个波源,R,i,（,i,=1,，,，,7,），,BC,边对等地安放,7,个接收器,S,j,（,j,=1,，,，,7,），记录由,R,i,发出的弹性波到达,S,j,的时间,ij,秒）,，,见表,2-3,。,已知弹性波在介质和空气中的传播速度分别为,2880,（米,/,秒）和,320,（米,/,秒），且弹性波沿板边缘的传播速度与在介质中的传播速度相同。,1,）确定该平板内空洞的位置。,2,）只根据由,P,i,发出的弹性波到达,Q,j,的时间,t,ij,（,i,，,j,=1,，,，,7,），能确定空洞的位置吗？讨论在同样能够确定空洞位置的前提下，减少波源和接受器的方法。,1.2,模型的假设及符号说明,（,1,）模型的假设波源和接收器的设置仅按图中的设置来考虑，不考虑其它情况。在图中的一个小方格要么全是空气，要么全。,是介质。,（,2,）,符号说明,n,+1,：,x,轴方向等距地放置波源,P,i,（,i,=17,）的个 数，相应地,n,是,x,轴方向的方格数（已知 的,n,=6,）；,m,+1：,y,轴方向等距地放置波源,R,i,（,i,=17）的个 数，相应地,m,是,y,轴方向的方格数（已知的 m=6）；,d,1,：,x,轴方向方格的边长（已知的，,d,1,=40米）；,d,2,：,y,轴方向方格的边长（已知的，,d,2,=40米）；,t,ij,（,i,=06；,j,=06）：由波源,P,i,发生的弹性波到达 接收器,Q,j,的时间；,p,ij,（,i=06,；,j=06,）：波线,P,i,Q,j,经历的介质的 总长度；,q,ij,（,i=06,；,j=06,）：波线,P,i,Q,j,经历的空洞的 总长度；,v,1,：弹性波在介质中的传播速度（已知的,v,1,=2880,米,/,秒）；,v,2,：弹性波在空气中的传播速度（已知的,v,2,=320,米,/,秒）；,1.3,问题的分析及数学模型,（,1,）问题的分析这是一个反问题，不妨先看它正问题。已知空洞的位置及弹性波在介质和空气中的传播速度，求由波源,P,i,发出的弹性波到达接收器,Q,j,的时间,t,ij,，为了求出,t,ij,，可先算出波线,P,i,Q,j,经历的每个小方格的长度，根据该方格是介质还是空洞，分别除以波在介质和空气中的传播速度，对经历的所有方格求和，即得,t,ij,。上述的正问题等价于：算出波线,P,i,Q,j,经历的各个小方格的长度，根据该方格是介质还是空洞，分别乘以,0,和,1,，对经历的所有方格求和，即可得到波线,P,i,Q,j,经历的空洞的长度。,反问题可以这样求解：先由波源,P,i,发出的弹性波到达接收器,Q,j,的时间,t,ij,及弹性波在介质和空气中的传播速度，算出波线,P,i,Q,j,经历的空洞的总长度，再由,P,i,Q,j,经历的每个方格的长度，即可解出每个方格是介质还是空洞（即是,0,还是,1,）。,（,2,）数学建模,建立坐标系如图,2-2,，,X,、,Y,轴上分别等距地放置,n,+1,，,m,+1,个波源（本题,n,=,m,=6,），将平板划分为 个方格，则,1,）波线,P,i,Q,j,经历的介质的总长度,p,ij,和空洞的总长度,q,ij,满足下面两个式子：,1,2,6,7,8,12,31,32,36,3,4,5,13,（如图22）,2,）波线,P,i,Q,j,的方程为,(3,）直线方程,(3),与 的交点可算出,从而求得,P,i,Q,j,经历的每个方格长度。,3,）将 个方格如图示顺序排列，第,i,个方格的特征量记为,(4,）全体方格的特征向量记为,。,4,）由,AB,至,CD,的波线,P,i,Q,j,共 条，其中,n,+1,条,P,i,Q,j,是无用的，去掉无用波线后,n(n,+1,）条波线,P,i,Q,j,的顺序按 ，；,；排列，每条波线经历的各个方格的长度排成一行向量，于是得到 条波线经历,mn,个方格的长度矩阵,5,）将,1,）得到的波线,P,i,Q,j,经历空洞的总长度,q,ij,，按照,4,）中波线的顺序构成空洞长度（列）向量 ，则应有 （,5,）,6,）完全类似地处理波线,R,k,S,l,，得到,m,（,m+,1,）条波线经历,mn,个方格的长度矩阵 和空洞长度（列）向量 ，构造,则,（6）,当,Rank(A)=mn 时可求出该方程最小二乘意义下的解,（7）,7,）取适当的 （,8,）取 的小方格为介质，小方格为空洞，其余的（如果有的话）无法确定。,2.3.4,模型求解,%,空洞探测,ch23.m,%,文件名：,ch23.m,%,计算,PiQj,在每个格子中的线段长度。,clear,n=1;,for i=0:6,for j=0:6,ai=i;aj=j;if i=j,aj=j+1;,if aj=7,break;,end,end,Chudian=ai,0;Modian=aj,6;Dianwei=ai,0;pp=;,for m=1:6,x=Chudian(1)+(Modian(1)-Chudian(1)*(m-Chudian(2)/.,(Modian(2)-Chudian(2);,Dianwei=Dianwei;x,m;end,if ajai,p=aj;aj=ai;ai=p;end,if abs(aj-ai)=2,for k=(aj+1):(ai-1),y=Chudian(2)+(Modian(2)-Chudian(2)*(k-Chudian(1)/.,(Modian(1)-Chudian(1);,pp=pp;k,y;end,end,Jiaodian=Dianwei;pp;hang,lie=size(Jiaodian);aa=Jiaodian(:,1);bb=unique(aa);hang1,lie1=size(bb);Zuobiao=;,for ii=1:hang1,for jj=1:hang,if Jiaodian(jj,1)=bb(ii)Zuobiao(ii,:)=Jiaodian(jj,:);break;end;end;end;,for pp=1:hang1-1,d=sqrt(sum(Zuobiao(pp,:)-Zuobiao(pp+1,:).2);,if j=2 ifdistance_PQ(n,:)=distance_PQ(n-1,:);,distance_PQ(n,:)=;,else,n=n+1;,end,else,n=n+1;,end,end,end,%,计算矩阵,A1,的秩,A1=distance_PQ;,disp(,矩阵,A1,的秩,),rank(A1),%,计算,RiSj,在每个格子中的线段长度。,n=1;,for i=0:6,for j=0:6,ai=i;aj=j;,if i=j,aj=j+1;,if aj=7,break;,end,end,Chudian=0,ai;Modian=6,aj;Dianwei=0,ai;pp=;,for m=1:6 y=Chudian(2)+(Modian(2)-Chudian(2)*(m-Chudian(1)/.,(Modian(1)-Chudian(1);,Dianwei=Dianwei;m,y;end,if ajai,p=aj;aj=ai;ai=p;end,if abs(aj-ai)=2,for k=(aj+1):(ai-1)x=Chudian(1)+(Modian(1)-Chudian(1)*(k-Chudian(2)/.,(Modian(2)-Chudian(2);,pp=pp;x,k;end,end,Jiaodian=Dianwei;pp;hang,lie=size(Jiaodian);aa=Jiaodian(:,2);bb=unique(aa);hang1,lie1=size(bb);,Zuobiao=;,for ii=1:hang1,for jj=1:hang,if Jiaodian(jj,2)=bb(ii),Zuobiao(ii,:)=Jiaodian(jj,:);break;end;end;end;,for pp=1:hang1-1,d=sqrt(sum(Zuobiao(pp,:)-Zuobiao(pp+1,:).2);,if j=2 ifdistance_RS(n,:)=distance_RS(n-1,:);,distance_RS(n,:)=;,else,n=n+1;,end,else,n=n+1;,end,end,end,%,对矩阵进行合并操作，计算系数矩阵,A,A=distance_PQ;distance_RS;,disp(,矩阵,A,的秩,),rank(A),%,首先输入时间矩阵,并给出,b,b1=0.0895 0.1996 0.2032 0.4181 0.4923 0.5646;,0.0989 0.4413 0.4318 0.4770 0.5242 0.3805;,0.3052 0.4132 0.4153 0.4156 0.3563 0.1919;,0.3221 0.4453 0.4040 0.1789 0.0740 0.2122;,0.3490 0.4529 0.2263 0.1917 0.1768 0.1810;,0.3807 0.3177 0.2364 0.3064 0.2217 0.1031;,0.4311 0.3397 0.3566 0.1954 0.0760 0.0688;,b2=0.0602 0.0813 0.3516 0.3867 0.4314 0.5721;,0.0753 0.2852 0.4341 0.3491 0.4800 0.4980;,0.3456 0.3205 0.4093 0.4240 0.4540 0.3112;,0.3655 0.3289 0.4247 0.3249 0.2134 0.1017;,0.3165 0.2409 0.3214 0.3256 0.1874 0.2130;,0.2749 0.3891 0.5895 0.3016 0.2058 0.0706;,0.4434 0.4919 0.3904 0.0786 0.0709 0.0914;,lie1=;lie2=;,for i=1:7,m=b1(i,:);n=b2(i,:);,lie1=lie1;m;lie2=lie2;n;,end,b=lie1;lie2;,%,计算方程,Ax=b,最小二乘意义下解,x=Ab;,tushi=;,for i=0:5 pp=x(i*6+1:(i+1)*6);tushi=tushi;pp;,end,x=tushi;,disp(,输出结果,X),x,%,进行空洞判断,a=0.11;b=-0.11;,for i=1:6,for j=1:6,if tushi(i,j)b&tushi(i,j)a,xx(i,j)=0;,else,xx(i,j)=1;,end,end,end,xx,执行后输出,矩阵,A1,的秩,ans=,29,矩阵,A,的秩,ans=,36,输出结果,X,x=,0.0098 0.0180 0.0058 0.0070 0.0181 0.0137,0.0145 0.1233 0.1325 0.0140 0.0137 0.0125,0.0203 0.1196 0.1204 0.0217 0.1155 0.0174,0.0255 0.0128 0.1273 0.1246 0.0173 0.0074,0.0127 0.1271 0.0145 0.0134 0.0069 0.0134,0.0054 0.0159 0.0134,0.0075 0.0185 0.0173,取,=0.1,有,xx=0 0 0 0 0,0,0 1 1 0 0 0,0 1 1 0 1 0,0 0 1 1 0 0,0 1 0 0 0 0,0 0 0 0 0 0,对应于平板上的空洞位置如图,23,（注意,:,相对于,x,的排列，上下正相反）。,图23,A,B,C,D,2.3.5,问题的进一步讨论,1,）只根据由,P,i,发出的弹性波到达,Q,j,的时间,t,ij,，可得长度矩阵 ，若,(5),式有解则能确定空洞的位置，但 ，,(5),式无解，故不能确定空洞的位置。,2,）在同样能够确定空洞位置的前提下，减少波源和接受器的方法：设在,BC,边仍安放,7,个接收器,S,l,，而,AD,边设置波源,R,k,的数量一个一个地增加（由下而上），发现设置两个波源时有 ，虽然（,7,）式有解，但由（,7,）式的解无法选取,得到 。直到设置,5,个波源时由（,7,）式才得到满意的解。,x,=-0.0324 0.0485 -0.0889 -0.0959 0.0231 0.0661 0.0181 0.9766 1.0591 0.0345 0.0201 -.0688 0.0563 0.9402 0.9846 0.1003 0.9100 -0.0468 0.0588 0.0430 1.0429 0.9694 0.0085 -0.0198 0.1416 0.9766 -0.1041 -0.0225 -0.0726 -.0070 -0.2188 0.0021 0.0664 -0.0444 0.0680 0.0376,取,=0.22,，即得到与前面同样的结果。另外，也可试验在,AB,和,AD,边的各,7,个点上，选取不同位置设置波源，并类似地考虑接收器，找出减少波源（和接收器）的方案。,2.1.3 MATLAB,求解,对于方程组 ，若,A,是满秩方阵，则,MATLAB,求解命令为：,MATLAB,中常用的矩阵运算命令,1,、求解方程组 若,A,是满秩的方阵，求解命令为：它等价于,2,、求解方程组 若,A,不是满秩的方阵，求解命令为：,rrefC,其中,C=A b3,、矩阵,A,的秩,rank,（,A,）,4,、,V,，,D=eig(A),输出的,D,是矩阵,A,的特征值构成的对角矩阵，输出的,V,是矩阵,A,的特征向量构成的特征向量矩阵；,