《数学建模非线性规划模型用MATLAB++LINGO课件.ppt》由会员分享，可在线阅读，更多相关《数学建模非线性规划模型用MATLAB++LINGO课件.ppt（68页珍藏版）》请在三一办公上搜索。
1、非线性规划,讲义大纲:,1.非线性规划的定义和相关概念.2.常用的求解非线性规划的方法.3.MATLAB求解非线性规划及例题.4.lingo求解非线性规划及例题.5.练习.,非现性规划的基本概念,定义 如果目标函数或约束条件中至少有一个是非线性函数时的最优化问题就叫做非线性规划问题,一般形式:（1）其中，是定义在 En 上的实值函数，简记:,其它情况:求目标函数的最大值或约束条件为小于等于零的情况，都可通过取其相反数化为上述一般形式,定义1 把满足问题（1）中条件的解 称为可行解（或可行点），所有可行点的集合称为可行集（或可行域）记为D即 问题(1)可简记为,定义2 对于问题(1)，设，若存在
2、,使得对一切，且，都有，则称X*是f(X)在D上的局部极小值点（局部最优解）特别地当 时，若则称X*是f(X)在D上的严格局部极小值点（严格局部最优解）,定义3 对于问题(1),设,对任意的,都有 则称X*是f(X)在D上的全局极小值点（全局最优解）特别地当 时，若，则称X*是f(X)在D上的严格全局极小值点（严格全局最优解）,返回,非线性规划的基本解法,SUTM外点法,SUTM内点法（障碍罚函数法）,1、罚函数法,2、近似规划法,返回,罚函数法,罚函数法基本思想是通过构造罚函数把约束问题转化为一系列无约束最优化问题，进而用无约束最优化方法去求解这类方法称为序列无约束最小化方法(Sequent
3、ial Unconstrained Minization Technique)简称为SUMT法其一为SUMT外点法，其二为SUMT内点法,其中T(X,M)称为罚函数，M称为罚因子，带M的项称为罚项，这里的罚函数只对不满足约束条件的点实行惩罚：当 时，满足各，故罚项=0，不受惩罚当 时，必有 的约束条件，故罚项0，要受惩罚,SUTM外点法,罚函数法的缺点是：每个近似最优解Xk往往不是容许解，而只能近似满足约束，在实际问题中这种结果可能不能使用；在解一系列无约束问题中，计算量太大，特别是随着Mk的增大，可能导致错误,1、任意给定初始点X0，取M11，给定允许误差，令k=1；2、求无约束极值问题 的
4、最优解，设为Xk=X(Mk)，即；3、若存在，使，则取MkM()令k=k+1返回（2），否则，停止迭代得最优解.计算时也可将收敛性判别准则 改为.,SUTM外点法(罚函数法)的迭代步骤,例题.min x2 s.t.x=1 显然本问题的最优解为x*=1,用SMT外点法:T(x,M)=x2+Mmin(0,x-1)2=求minT(x,M).本题可由 T(x,M)=2x+2M(x-1)=0,解得:x=M/(1+M),M趋于无穷.可知x从小于1趋于1,罚函数从外部趋于最优解.,SUTM内点法（障碍函数法）,内点法的迭代步骤,例题.min x2 s.t.x=1 显然本问题的最优解为x*=1,用SMT内点法
5、:罚函数I(x,rk)=x2 rkln(x-1),求其最值.可见,x从大于1趋于1.,近似规划法的基本思想：将问题(3)中的目标函数 和约束条件 近似为线性函数，并对变量的取值范围加以限制，从而得到一个近似线性规划问题，再用单纯形法求解之，把其符合原始条件的最优解作为(3)的解的近似,近似规划法,每得到一个近似解后，都从这点出发，重复以上步骤,这样，通过求解一系列线性规划问题，产生一个由线性规划最优解组成的序列，经验表明，这样的序列往往收敛于非线性规划问题的解。,近似规划法的算法步骤如下,返回,例题1、用近似规划法求解下列问题。,解：第一次迭代：在点x1处将g1(x),g2(x)线性化。,步长
6、限制：,解（1）（2）（3）（4），得，检验，可得x2不满足原始约束。,第二次迭代，减少,解（1,2,3,5)，得x2=(4,3)T,f(x2)=-11.所以x2为最小值点.,用MATLAB软件求解,其输入格式如下:1.x=quadprog(H,C,A,b);2.x=quadprog(H,C,A,b,Aeq,beq);3.x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB);4.x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB,X0);5.x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB,X0,options);6.x,fval=qua
7、prog(.);7.x,fval,exitflag=quaprog(.);8.x,fval,exitflag,output=quaprog(.);,1、二次规划,例1 min f(x1,x2)=-2x1-6x2+x12-2x1x2+2x22 s.t.x1+x22-x1+2x22 x10,x20,1、写成标准形式：,2、输入命令：H=1-1;-1 2;c=-2;-6;A=1 1;-1 2;b=2;2;Aeq=;beq=;VLB=0;0;VUB=;x,z=quadprog(H,c,A,b,Aeq,beq,VLB,VUB),3、运算结果为：x=0.6667 1.3333 z=-8.2222,s.t.
8、,1.首先建立M文件fun.m,定义目标函数F（X）:function f=fun(X);f=F(X);,2、一般非线性规划,其中X为n维变元向量，G(X)与Ceq(X)均为非线性函数组成的向量，其它变量的含义与线性规划、二次规划中相同.用Matlab求解上述问题，基本步骤分三步：,3.建立主程序.非线性规划求解的函数是fmincon,命令的基本格式如下：(1)x=fmincon(fun,X0,A,b)(2)x=fmincon(fun,X0,A,b,Aeq,beq)(3)x=fmincon(fun,X0,A,b,Aeq,beq,VLB,VUB)(4)x=fmincon(fun,X0,A,b,A
9、eq,beq,VLB,VUB,nonlcon)(5)x=fmincon(fun,X0,A,b,Aeq,beq,VLB,VUB,nonlcon,options)(6)x,fval=fmincon(.)(7)x,fval,exitflag=fmincon(.)(8)x,fval,exitflag,output=fmincon(.),输出极值点,M文件,迭代的初值,参数说明,变量上下限,注意：1 fmincon函数提供了大型优化算法和中型优化算法。默认时，若在fun函数中提供了梯度（options参数的GradObj设置为on），并且只有上下界存在或只有等式约束，fmincon函数将选择大型算法。当
10、既有等式约束又有梯度约束时，使用中型算法。2 fmincon函数的中型算法使用的是序列二次规划法。在每一步迭代中求解二次规划子问题，并用BFGS法更新拉格朗日Hessian矩阵。3 fmincon函数可能会给出局部最优解，这与初值X0的选取有关。,1、写成标准形式：s.t.,2x1+3x2 6 s.t x1+4x2 5 x1,x2 0,例2,2、先建立M-文件 fun3.m:function f=fun3(x);f=-x(1)-2*x(2)+(1/2)*x(1)2+(1/2)*x(2)2,3、再建立主程序youh2.m：x0=1;1;A=2 3;1 4;b=6;5;Aeq=;beq=;VLB=
11、0;0;VUB=;x,fval=fmincon(fun3,x0,A,b,Aeq,beq,VLB,VUB),4、运算结果为：x=0.7647 1.0588 fval=-2.0294,1先建立M文件 fun4.m,定义目标函数:function f=fun4(x);f=exp(x(1)*(4*x(1)2+2*x(2)2+4*x(1)*x(2)+2*x(2)+1);,例3,2再建立M文件mycon.m定义非线性约束：function g,ceq=mycon(x)g=x(1)+x(2);ceq=1.5+x(1)*x(2)-x(1)-x(2);-x(1)*x(2)-10;,3主程序youh3.m为:x0
12、=-1;1;A=;b=;Aeq=1 1;beq=0;vlb=;vub=;x,fval=fmincon(fun4,x0,A,b,Aeq,beq,vlb,vub,mycon),3.运算结果为：x=-1.2250 1.2250 fval=1.8951,例4,1先建立M-文件fun.m定义目标函数:function f=fun(x);f=-2*x(1)-x(2);,2再建立M文件mycon2.m定义非线性约束：function g,ceq=mycon2(x)g=x(1)2+x(2)2-25;x(1)2-x(2)2-7;,3.主程序fxx.m为:x0=3;2.5;VLB=0 0;VUB=5 10;x,f
13、val,exitflag,output=fmincon(fun,x0,VLB,VUB,mycon2),4.运算结果为:x=4.0000 3.0000fval=-11.0000exitflag=1output=iterations:4 funcCount:17 stepsize:1 algorithm:1x44 char firstorderopt:cgiterations:,应用实例：供应与选址,某公司有6个建筑工地要开工，每个工地的位置（用平面坐标系a，b表示，距离单位：千米）及水泥日用量d(吨)由下表给出。目前有两个临时料场位于A(5,1)，B(2,7)，日储量各有20吨。假设从料场到工地
14、之间均有直线道路相连。（1）试制定每天的供应计划，即从A，B两料场分别向各工地运送多少吨水泥，使总的吨千米数最小。（2）为了进一步减少吨千米数，打算舍弃两个临时料场，改建两个新的，日储量各为20吨，问应建在何处，节省的吨千米数有多大？,（一）、建立模型,记工地的位置为(ai，bi)，水泥日用量为di，i=1,6;料场位置为(xj，yj)，日储量为ej，j=1,2；从料场j向工地i的运送量为Xij。,当用临时料场时决策变量为：Xij，当不用临时料场时决策变量为：Xij，xj，yj。,（二）使用临时料场的情形,使用两个临时料场A(5,1)，B(2,7).求从料场j向工地i的运送量为Xij，在各工地
15、用量必须满足和各料场运送量不超过日储量的条件下，使总的吨千米数最小，这是线性规划问题.线性规划模型为：,设X11=X1,X21=X 2,X31=X 3,X41=X 4,X51=X 5,X61=X 6X12=X 7,X22=X 8,X32=X 9,X42=X 10,X52=X 11,X62=X 12 编写程序gying1.m,gying1.m,clear a=1.25 8.75 0.5 5.75 3 7.25;b=1.25 0.75 4.75 5 6.5 7.75;d=3 5 4 7 6 11;x=5 2;y=1 7;e=20 20;for i=1:6 for j=1:2 aa(i,j)=sqr
16、t(x(j)-a(i)2+(y(j)-b(i)2);endendCC=aa(:,1);aa(:,2);,A=1 1 1 1 1 1 0 0 0 0 0 0;0 0 0 0 0 0 1 1 1 1 1 1;B=20;20;Aeq=1 0 0 0 0 0 1 0 0 0 0 0%从第一二料场运到工地一的料 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 1;beq=d(1);d(2);d(3);d(4);d(5)
17、;d(6);VLB=0 0 0 0 0 0 0 0 0 0 0 0;VUB=;x0=1 2 3 0 1 0 0 1 0 1 0 1;xx,fval=linprog(CC,A,B,Aeq,beq,VLB,VUB,x0),计算结果为：,x=3.0000 5.0000 0.0000 7.0000 0.0000 1.0000 0.0000 0.0000 4.0000 0.0000 6.0000 10.0000fval=136.2275,（三）改建两个新料场的情形,改建两个新料场，要同时确定料场的位置(xj,yj)和运送量Xij，在同样条件下使总吨千米数最小。这是非线性规划问题。非线性规划模型为：用li
18、ngo解此题 MODEL:,设 X11=X1,X21=X 2,X31=X 3,X41=X 4,X51=X 5,X61=X 6 X12=X 7,X22=X 8,X32=X 9,X42=X 10,X52=X 11,X62=X 12 x1=X13,y1=X14,x2=X15,y2=X16,（1）先编写M文件liaoch.m定义目标函数。,(2)取初值为线性规划的计算结果及临时料场的坐标:x0=3 5 0 7 0 1 0 0 4 0 6 10 5 1 2 7;编写主程序gying2.m.,function f=liaoch(x)a=1.25 8.75 0.5 5.75 3 7.25;b=1.25 0.
19、75 4.75 5 6.5 7.75;d=3 5 4 7 6 11;e=20 20;f1=0;for i=1:6 s(i)=sqrt(x(13)-a(i)2+(x(14)-b(i)2);f1=s(i)*x(i)+f1;endf2=0;for i=7:12 s(i)=sqrt(x(15)-a(i-6)2+(x(16)-b(i-6)2);f2=s(i)*x(i)+f2;endf=f1+f2;,clear%x0=3 5 0 7 0 1 0 0 4 0 6 10 5 1 2 7;%x0=3.0000 5.0000 0.0707 7.0000 0 0.9293 0 0 3.9293 0 6.0000 1
20、0.0707 6.3875 4.3943 5.7511 7.1867;%x0=3.0000 5.0000 0.3094 7.0000 0.0108 0.6798 0 0 3.6906 0 5.9892 10.3202 5.5369 4.9194 5.8291 7.2852;x0=3 5 4 7 1 0 0 0 0 0 5 11 5.6348 4.8687 7.2479 7.7499;A=1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 0 0 0 0;B=20;20;,原程序gying2.m,Aeq=1 0 0 0 0 0 1 0
21、0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0;beq=3 5 4 7 6 11;vlb=zeros(12,1);-inf;-inf;-inf;-inf;vub=;x,fval,exitflag=fmincon(liaoch,x0,A,B,Aeq,beq,vlb,vub),(3)计算结果为：,x
22、=3.0000 5.0000 0.0707 7.0000 0 0.9293 0 0 3.9293 0 6.0000 10.0707 6.3875 4.3943 5.7511 7.1867fval=105.4626exitflag=1,(4)若修改主程序gying2.m,取初值为上面的计算结果:x0=3.0000 5.0000 0.0707 7.0000 0 0.9293 0 0 3.9293 0 6.0000 10.0707 6.3875 4.3943 5.7511 7.1867,得结果为:x=3.0000 5.0000 0.3094 7.0000 0.0108 0.6798 0 0 3.69
23、06 0 5.9892 10.3202 5.5369 4.9194 5.8291 7.2852fval=103.4760exitflag=1,总的吨千米数比上面结果略优.,(5)若再取刚得出的结果为初值,却计算不出最优解.,(6)若取初值为:x0=3 5 4 7 1 0 0 0 0 0 5 11 5.6348 4.8687 7.2479 7.7499,则计算结果为:x=3.0000 5.0000 4.0000 7.0000 1.0000 0 0 0 0 0 5.0000 11.0000 5.6959 4.9285 7.2500 7.7500fval=89.8835exitflag=1总的吨千米
24、数89.8835比上面结果更好.gying2.m,通过此例可看出fmincon函数在选取初值上的重要性.用lingo解此题的方法:,Lingo的使用数学函数,LINGO提供了大量的标准数学函数：abs(x)返回x的绝对值sin(x)返回x的正弦值，x采用弧度制cos(x)返回x的余弦值tan(x)返回x的正切值exp(x)返回常数e的x次方log(x)返回x的自然对数lgm(x)返回x的gamma函数的自然对数sign(x)如果x=0时，返回不超过x的最大整数；当x0时，返回不低于x的最大整数。smax(x1,x2,xn)返回x1，x2，xn中的最大值smin(x1,x2,xn)返回x1，x2
25、，xn中的最小值,例4.3 给定一个直角三角形，求包含该三角形的最小正方形。,LINGO代码如下：model:Title trangle;sets:object/1.3/:f;endsetsdata:a,b=3,4;!两个直角边长，修改很方便;Enddataf(1)=a*sin(x);f(2)=b*cos(x);f(3)=a*cos(x)+b*sin(x);min=smax(f(1),f(2),f(3);bnd(0,x,1.57);end,例题:,一奶制品加工厂用牛奶生产A1,A2两种奶制品，1桶牛奶可以在甲车间用12小时加工成3公斤A1，或者在乙车间用8小时加工成4公斤A2。根据市场需求，生
26、产的A1,A2全部能售出，且每公斤A1获利24元，每公斤A2获利16元。现在加工厂每天能得到50桶牛奶的供应，每天正式工人总的劳动时间480小时，并且甲车间每天至多能加工100公斤A1，乙车间的加工能力没有限制。试为该厂制订一个生产计划，使每天获利最大，并进一步讨论以下3个附加问题：1）若用35元可以买到1桶牛奶，应否作这项投资？若投资，每天最多购买多少桶牛奶？2）若可以聘用临时工人以增加劳动时间，付给临时工人的工资最多是每小时几元？3）由于市场需求变化，每公斤A1的获利增加到30元，应否改变生产计划？,模型代码如下：max=72*x1+64*x2;x1+x2=50;12*x1+8*x2=48
27、0;3*x1=100;,求解这个模型,结果如下:Global optimal solution found at iteration:0 Objective value:3360.000 Variable Value Reduced Cost X1 20.00000 0.000000 X2 30.00000 0.000000 Row Slack or Surplus Dual Price 1 3360.000 1.000000 2 0.000000 48.00000 3 0.000000 2.000000 4 40.00000 0.000000,Ranges in which the basi
28、s is unchanged:Objective Coefficient Ranges Current Allowable Allowable Variable Coefficient Increase Decrease X1 72.00000 24.00000 8.000000 X2 64.00000 8.000000 16.00000 Righthand Side Ranges Row Current Allowable Allowable RHS Increase Decrease 2 50.00000 10.00000 6.666667 3 480.0000 53.33333 80.0
29、0000 4 100.0000 INFINITY 40.00000,做灵敏性分析，结果如下:,用lingo解此题MODEL:Title Location Problem;sets:demand/1.6/:a,b,d;supply/1.2/:x,y,e;link(demand,supply):c;endsetsdata:!locations for the demand(需求点的位置);a=1.25,8.75,0.5,5.75,3,7.25;b=1.25,0.75,4.75,5,6.5,7.75;!quantities of the demand and supply（供需量）;d=3,5,4,
30、7,6,11;e=20,20;enddata,init:!initial locations for the supply（初始点）;x,y=5,1,2,7;endinit!Objective function（目标）;OBJ min=sum(link(i,j):c(i,j)*(x(j)-a(i)2+(y(j)-b(i)2)(1/2);!demand constraints（需求约束）;for(demand(i):DEMAND_CON sum(supply(j):c(i,j)=d(i););!supply constraints（供应约束）;for(supply(i):SUPPLY_CON s
31、um(demand(j):c(j,i)=e(i););for(supply:bnd(0.5,X,8.75);bnd(0.75,Y,7.75););END,运行的结果为:,X(1)5.695967X(2)7.249997 Y(1)4.928560 Y(2)7.750000Feasible solution found.Objective value:89.88349 Total solver iterations:81,1.3 应用举例(二)-整数规化,解 先看有多少种裁料方案，再进行组合和选择。方案：,例1 合理利用线材问题 现要做一百套钢管,每套要长为2.9m、2.1m和1.5m的钢管各一根
32、。已知原料长7.4m，问应如何下料，使用的原料最省。,设用方案,分别裁原料钢管x1,x2,x8根,则：,上题中,如果每套1.8m的钢管要70根,要求使用的切割模式不超过3种.问应如何下料，使用的原料最省。,解:设xi按第i种模式切割的原料钢管根数(i=1,2,3)r1i,r2i,r3i,r4i第i种切割模式下,每根原料生产2.9m、2.1m和1.5m,1.8m的钢管的数量.目标函数(总根数)(套数约束),按第i种切割模式下,每根钢管的长度限制:,这是一个整数非线性规划模型，用LINGO要运行很长时间也难以得到最优解。,因三种切割模式的排列顺序无关紧要，所以不妨增加以下约束：,又钢管的总根数有明
33、显的上界和下界。首先，原料钢管的总根数不可能少于,其次，考虑第一种切割模式下只生产2.9m钢管，一根原料钢管切割成2根2.9m钢管，100套钢管需要50根原料；如此可计算出钢管的上界：,所以可以增加以下约束：,Lingo源程序,model:title 钢管下料；min=x1+x2+x3;x1*r11+x2*r12+x3*r13=100;x1*r21+x2*r22+x3*r23=100;x1*r31+x2*r32+x3*r33=100;x1*r41+x2*r42+x3*r43=70;x1+x2+x3=105;x1+x2+x3=152;,2.9*r11+2.1*r21+1.5*r31+1.8*r4
34、1=5.9;2.9*r12+2.1*r22+1.5*r32+1.8*r42=5.9;2.9*r13+2.1*r23+1.5*r33+1.8*r43=5.9;2.9*r11+2.1*r21+1.5*r31+1.8*r41=x2;x2=x3;gin(x1);gin(x2);gin(x3);gin(r11);gin(r21);gin(r31);gin(r41);gin(r12);gin(r22);gin(r32);gin(r42);gin(r13);gin(r23);gin(r33);gin(r43);end,运行结果如下：,Feasible solution found at iteration:
35、26565 Model Title:钢管下料；Variable Value X1 72.00000 R11 1.000000 X2 29.00000 R12 0.000000 X3 14.00000 R13 2.000000 R21 1.000000 R22 1.000000 R23 0.000000 R31 0.000000 R32 3.000000 R33 1.000000 R41 1.000000 R42 0.000000 R43 0.000000,Lingo源程序二,model:Title 钢管下料-最小化钢管根数的LINGO模型;SETS:NEEDS/1.4/:LENGTH,NUM;
36、!定义基本集合NEEDS及其属性LENGTH,NUM;CUTS/1.3/:X;!定义基本集合CUTS及其属性X;PATTERNS(NEEDS,CUTS):R;!定义派生集合PATTERNS（这是一个稠密集合）及其属性R;ENDSETSDATA:LENGTH=2.9 2.1 1.5 1.8;NUM=100 100 100 70;CAPACITY=7.4;ENDDATAmin=SUM(CUTS(I):X(I);!目标函数;,Lingo源程序二,FOR(NEEDS(I):SUM(CUTS(J):X(J)*R(I,J)NUM(I);!满足需求约束;FOR(CUTS(J):SUM(NEEDS(I):LE
37、NGTH(I)*R(I,J)CAPACITY-MIN(NEEDS(I):LENGTH(I);!合理切割模式约束;SUM(CUTS(I):X(I)105;SUM(CUTS(I):X(I)X(I+1);!人为增加约束;FOR(CUTS(J):GIN(X(J);FOR(PATTERNS(I,J):GIN(R(I,J);end,运行结果如下：,Local optimal solution found at iteration:64408 Objective value:113.0000 Model Title:钢管下料-最小化钢管根数的LINGO模型 Variable Value Reduced Co
38、st CAPACITY 7.400000 0.000000 LENGTH(1)2.900000 0.000000 LENGTH(2)2.100000 0.000000 LENGTH(3)1.500000 0.000000 LENGTH(4)1.800000 0.000000 NUM(1)100.0000 0.000000 NUM(2)100.0000 0.000000 NUM(3)100.0000 0.000000 NUM(4)70.00000 0.000000 X(1)50.00000 1.000000 X(2)37.00000 1.000000 X(3)26.00000 1.000000,
39、运行结果如下：,R(1,1)1.000000 0.000000 R(1,2)0.000000 0.000000 R(1,3)2.000000 0.000000 R(2,1)2.000000 0.000000 R(2,2)0.000000 0.000000 R(2,3)0.000000 0.000000 R(3,1)0.000000 0.000000 R(3,2)2.000000 0.000000 R(3,3)1.000000 0.000000 R(4,1)0.000000 0.000000 R(4,2)2.000000 0.000000 R(4,3)0.000000 0.000000,某厂向用户提供发动机，合同规定，第一、二、三季度末分别交货40台、60台、80台每季度的生产费用为（元），其中x是该季生产的台数若交货后有剩余，可用于下季度交货，但需支付存储费，每台每季度c元已知工厂每季度最大生产能力为100台，第一季度开始时无存货，设a=50、b=0.2、c=4，问工厂应如何安排生产计划，才能既满足合同又使总费用最低讨论a、b、c变化对计划的影响，并作出合理的解释,练习,解：设每季度生产的发动机数为，则目标函数为,约束条件为：,
