1、Matlab优化,第 三 讲,本章主要学习利用matlab解决优化问题: 一、函数和方程1、函数(定义与求值)2、多项式(运算与求根)3、解方程(组) 二、优化问题1、无约束优化(函数的极值)2、规划问题(线性和非线性)3、目标规划4、数学模型,一、函数和方程,1、函数(1) 函数定义方法Fun=inline(funstr,var)%定义一个lnline函数,funstr是函数的表达式,var是变量名.Fun=(var)funstr%定义一个匿名函数,funstr是函数的表达式,var是变量名.Fun=Mfun%定义一个函数句柄,Mfun是函数的m文件名.,(2) 函数值的求法按上述方法定义好
2、函数后,求函数值就有下述两种方法:Fun(value)%求Fun在value处的值,value可以是一个数、向量或矩阵.feval(Fun,value)%feval是matlab内部命令,求Fun在value处的值,value可以是一个数、向量或矩阵.例1 求y=x*sin(x)在x=pi/4,1,pi/3,3/2处的值.clearfun1=inline(x.*sin(x),x);,fun2=(x) x.*sin(x); fun3=mfun其中 mfun.m 的语句如下:function y=mfun(x)y= x.*sin(x); x=pi/4,1,pi/3,3/2 y1=fun1(x);
3、y2=fun2(x); y3=fun3(x); z1=feval(fun1,x); z2=feval(fun2,x);z3=feval(fun3,x); 求出结果完全一样,为 0.55536 0.84147 0.9069 1.4962,2、多项式(1) 定义 matlab中一个多项式用系数降幂排列向量来表示.如:多项式x3+2x2-5,在matlab中表示为1 2 0 -5 (注意,不要遗漏系数为0的项).(2) 运算: 1) 求函数值 y=polyval(p,x)%求多项式p在x处的值.2) 乘法 p=conv(p1,p2)%p返回p1与p2的乘积.3) 除法 p,r=deconv(p1,p
4、2)%p返回p1除以p2的商,r返回余项.,(3) 求根 x=roots(p)% 求多项式p的根(复数范围内求解). p1=1 2 3 0 4 1; p2=1 0 1; p3=conv(p1,p2); p4,r=deconv(p1,p2); y=ployval(r,1); x=roots(r); 求解结果为:p3 = 1 2 4 2 7 1 4 1p4 = 1 2 2 -2r = 2 3x = -1.5000,3、解方程(组),matlab指令(下面两个命令只能求实根)x=fzero(Fun,x0)%返回一元函数Fun的一个零点,其中Fun是按前面方式定义的函数. X0为标量时,返回函数在x0
5、附近的零点;x0为区间a,b时,fzero 要求Fun在区间端点的函数值异号(否则fzero无法求解),返回函数在区间a,b中的零点.x,f,h=fsolve(Fun,x0)%返回一元或多元函数Fun的一个零点,其中Fun是按前面方式定义的函数. X0为迭代初值,f 返回函数Fun在x的函数值, 应该接近0;h返回值如果大于0, 说明计算结果可靠, 否则不可靠.,(2) 解方程(组),例1 求函数y=xsin(x2-x-1)在(-2,-0.1)内的零点.分析 首先要定义函数;然后考虑fzero是否可以 求解(因为给定的是一个区间),若不可求解,则用 fsolve求解. fun=inline(x
6、*sin(x2-x-1),x); fun(-2);fun(-0.1);%端点函数值是否异号可以看出端点函数直同号,fzero不可直接求解. 此时,我们作出函数图形,观察出一个初值来. fplot(fun,-2,-0.1); %观察图形得初值 x1=fsolve(fun,x01); x2=fsolve(fun,x02);,从右图可以看出 函数在x01=-1.6 及x02=-0.6附近 各有一个零点, 所以以此为初值 用fzero求根.,求解结果:x1=-1.5956 (x01=-1.6为初值)x2=-0.6180 (x02=-0.6为初值) 用fsolve求得相同的解. 从上述求解过程看出,有时
7、候需要“数形”结合才能求得解,这是由matlab求解方程(组)对初值的高度依赖性决定的.,例2 求解方程组,求方程组在原点附近的解:分析 首先要定义函数;然后用fzero或fsolve求解. 此处我们用inline定义函数 fun=inline(4*x(1)-x(2)+exp(x(1)/10-1,-x(1)+4*x(2)+x(1)2/8,x); x,f,h=fsolve(fun,0 0) 求解得:x = 0.2326 0.0565f = 1.0e-006 * 0.0908 0.1798h = 1 %解可靠,二、优化问题,1、无约束优化函数的极值除了自变量限于定义域外没有别的限制条 件的优化称为
8、无约束优化,典型的无约束优化 就是求函数的极值.min(y) 返回向量y的最小值;max(y) 返回向量y的最大值;x,f=fminbnd(fun,a,b) % x 返回一元函 数y=f(x)在a,b内局部极小值点,f 返回局部 极小值. 其中 fun 由前定义.,x,f=fminsearch(fun,x0) % x 返回一元 或多元函数 y=f(x) 在初值x0附近的局部极小 值点, f 返回局部极小值. 这里x,x0均为向量.x,f=fminunc(fun,x0) % x 返回一元 或多元函数 y=f(x) 在初值x0附近的局部极小 值点, f 返回局部极小值. 这里x,x0均为向量. 注
9、 当函数的阶数大于2时,使用 fminunc 比 fminsearch更有效,但当函数高度不连续时, 使用fminsearch更有效. 例1 求y=x*sin(x2-x-1)在-2,0上的极小值.,求解语句如下: (1) 用 min 命令. fun=inline(x*sin(x2-x-1),x) x=-2:0.01:0; m,k=min(fun(x);%返回最小值及其编址. 精度不高! x(k) %返回最小值点. 精度不高! (2) 用 fminbnd 命令. x,f=fminbnd(fun,-2,0) 结果为:x =-1.2455;f =-1.2138 (3) 用 fminsearch 命令
10、.,先“图形”结合,观察出初值x0=-1,然后求解. x,f=fminsearch(fun,-1) 结果为:x =-1.2455;f =-1.2138. 注 fminsearch 对初值x0有高度的依赖性! (4) 用 fminunc 命令. x,f=fminunc(fun,-1) 结果为:x =-1.2455;f =-1.2138. 注 fminunc 也对初值x0有高度的依赖性!例2 求z=5-x4-y4+4xy在原点附近极大值.,fun=inline(x(1)4+x(2)4-4*x(1)*x(2)-5,x) x,f=fminsearch(fun,0 0) 目标函数值(最大): 7 x:
11、0.999999999795223 -1.00000000022399例3 求z=sin(sqrt(x2+y2)/sqrt(x2+y2)在原 点附近的极大值. -pi/2fun=inline(sin(sqrt(x(1)2+x(2)2)/sqrt(x(1)2+x(2)2),x); x,f,h=fminsearch(fun,0 0)求出的结果为: x = 4.4180, 0.8198; f =-0.2172; h=1.可以看出,求出的误差很大,因为很容易看出, 在原点处取得最大值1.,自变量除了受到定义域的限制外,还有其他 条件的约束的优化问题称为规划问题. 本目介 绍matlab求解规划问题.
12、(1) 线性规划定义 目标函数和约束条件全为线性的规 划,其标准格式为:min fx s.t. Ax=b 其中 f 为目标函数的系数列向量,A为约束条 件的系数矩阵,b 为约束条件的右端项.,2、规划问题线性、二次、非线性,(LP0),还有一些衍生形式:min fxs.t. Ax=bAeq x=beq %等式约束min fx s.t. Ax=bAeq x=beqlb=x=ub %自变量范围min fxs.t. Ax=bAeq x=beqlb=x=ubx0 %迭代初值,(LP1),(LP2),(LP3),matlab 指令 x=linprog(f,A,b) %求解线性规划(LP0); x=lin
13、prog(f,A,b,Aeq,beq) %求解线性规划(LP1); x=linprog(f,A,b,Aeq,beq,lb,ub) %求解线性规划(LP2); x=linprog(f,A,b,Aeq,beq,lb,ub,x0) %求解线性规划(LP3);在上述各指令中,可以用代替缺省参数;左端 输出项若为x,Fval,则可得最优解处的最优值.,例1 解下列线性规划问题:,f=-5;-4;-6; %目标函数的系数列向量 A=1 -1 1;3 2 4;3 2 0; %约束条件的系数矩阵 b=20;42;30; %约束条件的右端项 lb=zeros(3,1); %自变量取值的下界 x,fval=lin
14、prog(f,A,b,lb,) 求解结果:x=0.0000 15.0000 3.0000; fval=-78.0000.,(2) 二次规划,定义 目标函数为二次多项式、约束条件全为 线性的规划问题. 其标准形式为:min xHx/2 + fxs.t. Ax=bAeq x=beqlb=x=ub 其中H,A为矩阵, f,x,b,lb,ub为列向量. matlab 指令 x,fval=quadprog(H,f,A,b,Aeq,beq,lb,ub,x0) 其中x0为迭代初值.,例2 求解下面二次规划:,H=1,-1;-1,2; %目标函数中二次型的矩阵 f=-2;-6; %目标函数中线性型的系数列向量
15、 A=1 1;-1 2;2 1; %约束条件的系数矩阵 b=2;2;3; %约束条件的右端项 lb=zeros(2,1); %自变量取值的下界 x,fval=quadprog(H,f,A,b,lb,) 求解结果:x=0.6667 1.3333; fval= -8.2222.,(3) 非线性规划,1) 标准型min f(x) %目标函数 s.t. Ax=b %线性不等约束Aeq x=beq %线性等式约束C(x)=0 %非线性不等约束Ceq(x)=0 %非线性等式约束 lb=x=ub %自变量的上下界 2) matlab 指令x,fval=fmincon(fun,初值,A,b,Aeq,beq,l
16、b,ub,nonlincon),对命令中的各参数解释一下: fun目标函数,一般用m函数文件定义; nonlincon非线性约束(包括等式约束和不等式约束),一般也用m函数文件定义. 例3 解下列非线性规划问题:,解 首先编写目标函数fun和非线性约束函数,nonlincon如下. 目标函数fun function f=fun(x) f=-x(1)2*x(2)*x(3)2/(2*x(1)3*x(3)2+3*x(1)2*x(2)2+2*x(2)2*x(3)3+x(1)3*x(2)2*x(3)2); 非线性约束函数nonlincon function c,ceq=nonlincon(x) %非线性
17、不等约束 c=1-x(1)2-x(2)2-x(3)2;x(1)2+x(2)2+x(3)2-4;,然后调用命令fmincon求解: x,fval=fmincon(fun,1 1 1,0 0 0,nonlincon) 结果为:x = 0.8676 0.5330 1.3107fval = -0.1537.,在实际中,人们碰到的大部分问题都是求某个目标函数的最大值或最小值的问题。但在某些情形下,人们会碰到一类特殊的问题,这些问题要求使最大值最小化才有意义。例如,某城市规划中需要确定急救中心、消防中心等的位置,可取的目标函数应该是到所有地点最大距离的最小值,而不是所有目的地的距离的和为最小。这是两种完全
18、不同的准则,在实际中也有广泛的应用。称之为最大最小化问题。1) 定义 极小化极大(minmax)问题的标准形式为:,其中x,lb,ub是向量,F(x)=Fi (x)和G(x)是函数向量.2) matlab 指令 x,fval,maxfval=fminimax(fun,x0,A,b,Aeq,beq,lb,ub,nonlincon) 参数说明:右端各参数完全同于 fmincon; x-最优点;fval-最优点处的值;maxfval-目标函数在x处的最大值.,例4 求下列函数最大值的最小化问题,解 先建立目标函数,并保存为myfun.m.然后在命令窗口输入如下命令: x0=0.1; 0.1 %初值,
19、x,fval,maxfval=fminimax(myfun,x0) 结果为:x=4.0000 4.0000; mavfval=0; fval=0 -64.0000 -2.0000 -8.0000 -0.0000. myfun.m附后: function f,g=myfun(x) f(1)=2*x(1)2+x(2)2-48*x(1)-40*x(2)+304; f(2)=-x(1)2-3*x(2)2; f(3)= x(1)+3*x(2)-18; f(4)=-x(1)-x(2); f(5)= x(1)+x(2)-8; g=; %无约束,例5 求例4中各函数的绝对值的最大值最小化问题.,解 应用例4的
20、myfun.m, 然后运行下述命令: x0=0.1 0.1; options=optimset(MinAbsMax,5); x,fval,maxfval=fminimax(myfun,x0,options) 结果为: x=4.9256 2.0796 maxfval=37.2356 fval=37.2356 -37.2356 -6.8357 -7.0052-0.9948,3、目标规划,线性规划的缺点: 1、目标单一性;2、约 束条件相容性;3、约束条件“刚性”.目标规划就是在线性规划的基础上为克服 线性规划的缺点和适应复杂的多目标最优决 策的需要发展起来的,它对众多的目标分别 确定一个希望实现的
21、目标值,然后按目标的 重要级别依次进行考虑和计算,以求得最接 近实现各目标预定数值的方案。目标规划特 别适合于经济活动中的目标管理。,matlab求解目标规划: 1、matlab中目标规划的标准型:min 目标函数(多个)s.t. 每个目标函数=goalc(x)=0 %非线性不等约束ceq(x)=0 %非线性等式约束Ax=b %线性不等约束Aeq x=beq %线性等式约束lb=x=ub 2、命令:fgoalattain格式:多种调用格式,最详细的调用格式为:x, fval, attainfactor, exitflag= fgoalattain( fun, x0, goal, weight,
22、 A, b,Aeq, beq, lb, ub, nonlincon, options ) 其中, attainfactor 为目标达到因子:,exitflag的值描述计算退出原因,此处忽略.,应用举例:,例6 生产安排问题 某企业生产甲乙两种产品,需要用到ABC三 种设备,关于产品的盈利与使用设备的工时 及限制如下表。,问:该企业应如何安排生产,使得在计划其内 (1) 力求使利润指标不低于1500元; (2) 考虑到市场需求,甲乙两种产品的产量尽量 比保持1:2; (3) 设备A为贵重设备,严格禁止超时使用; (4) 设备C可以适当加班,但要控制;设备B既要求充分利用,又尽可能不加班;在重要性
23、上,B是C的3倍.解 此问题明显非线性规划所能解决. 设生产,甲x1件,乙x2件。则由题意,有 (1) 利润 200x1+300x2+dm1-dp1=1500; (2) 产量比例 2x1-x2+dm2-dp2=0; (3) 设备A严禁超时 2x1+2x2=12; (4) 设备B 4x1+dm3-dp3=16; (5) 设备C 5x2+dm4-dp4=15; (6) 考虑设备B、C的重要性 目标表示为: (1) min dm1 (利润)(2) min dm2+dp2 (产量)(3) min dm3+dp3 (B)(4) min dp4 (C)(5) min 3(dm3+dp3)+dp4 (重要性
24、),应用matlab求解,建立了2个模型,见程序链接.需要说明的是,求解目标规划lingo与matlab 各有千秋. 同样的问题, 利用lingo建模得到的 结果与matlab有明显不同,而且lingo结果要好 一些, 见程序运行结果. 引起误差的主要原因 是各目标的优先级处理的不同.关于整数规划,目前matlab只能处理0-1线 性规划, 因为lindo或lingo在线性规划方面的 优势,在此我们就不介绍matlab求解了.,4、数学模型,本目介绍关于优化问题的数学模型及其求解. (1) 购房贷款的利率住房是居民消费的一个重要部分. 大部分 人选择银行按揭贷款,然后在若干年内逐月分 期还款.
25、 如果你借了10万,还款额一定超过10 万. 下面我们对购房贷款作一点分析.设贷款总额为x0, 贷款期限为n个月, 采取 逐月等额方式偿还本息. 若 x(k) 为第k个月的 欠款数, a为月还款数,r为月利率. 我们得到下,列递推关系式 X(k+1)=(1+r)x(k)-a, 那么, 递推可得如下公式: x(k)=(1+r)kx0-a(1+r)k-1/r; 由此得到月还款额度为: a=(1+r)nrx0/(1+r)n-1. 模型应用 下面是某报2000年3月30日的一 则房产广告:,不难算出,你向银行总共贷了25.2万,30年内,(*),共要还51.696万,是贷款额的2倍多. 这个案例 中贷
26、款年利率是多少呢?有人或许会这样计算: 年利率=(51.696-25.2)/30/25.2=3.5%, 但这个结果是错误的,因为你并不是等到30年 后一次性还款.根据a=0.1436,x0=25.2,n=360,由模型(*) 得到关于月利率r的高次代数方程: 25.2r(1+r)360-0.1436(1+r)360-1=0. 理论上可以用roots求解,但现实上并不合适.,我们用fzero求解. 以活期月利率为迭代初值, 即 r0=0.0198/12. 还要剔除r=0这个无意义 的解,在求解时我们将方程变了形: clear r0=0.0198/12; lilv=inline(25.2*(1+r
27、)360- 0.1436*(1+r)360-1)/r,r); r=fzero(lilv,r0); %月利率 R=12*r %年利率 求解结果为:r=0.0046; R=0.0553 (5.53%).,(2) 横渡长江,2002年5月1日, 横渡长江竞赛在武汉举行.假设在竞赛区域两岸为平行直线,它们之间 的垂直距离为 H=1160m, 从武昌汉阳门的正 对岸到汉阳南岸咀的距离为 L=1000m, 汉江 的平均流速为1.89m/s. 见下图:,1) 试为一个速度保持在1.5m/s的人选择游泳 方向,并估计他的成绩; 2) 若流速沿离岸边距离的分布为(社从武昌 汉阳门垂直向上为y轴正向):且游泳者的
28、速度大小仍全程保持不变,为1.5m /s,试为他选择游泳方向和路线,估计他的成绩. 3) 若流速沿离岸边距离为连续分布,例如,你又如何处理上述问题?,小结:matlab优化命令集锦,Optimization Toolbox - Version 3.1.1 (R2007a) 01-Feb-2007 1、Nonlinear minimization of functions.fminbnd - Scalar bounded nonlinear function minimization.fmincon - Multidimensional constrained nonlinear minimiz
29、ation.fminsearch - Multidimensional unconstrained nonlinear minimization, by Nelder-Mead direct search method.fminunc - Multidimensional unconstrained nonlinear minimization.fseminf - Multidimensional constrained minimization, semi-infinite constraints.,2、Nonlinear minimization of multi-objective fu
30、nctions.fgoalattain - Multidimensional goal attainment optimization fminimax - Multidimensional minimax optimization.3、Linear least squares (of matrix problems).lsqlin - Linear least squares with linear constraints.lsqnonneg - Linear least squares with nonnegativity constraints.4、Nonlinear least squ
31、ares (of functions).lsqcurvefit - Nonlinear curvefitting via least squares (with bounds).lsqnonlin - Nonlinear least squares with upper and lower bounds.,5、Nonlinear zero finding (equation solving).fzero - Scalar nonlinear zero finding.fsolve - Nonlinear system of equations solve (function solve).6、
32、Minimization of matrix problems.bintprog - Binary integer (linear) programming.linprog - Linear programming.quadprog - Quadratic programming.7、Controlling defaults and options.optimset - Create or alter optimization OPTIONS structure. optimget - Get optimization parameters from OPTIONS structure.,8、
33、Graphical user interface and plot routinesoptimtool - Optimization Toolbox Graphical User Interfaceoptimplotconstrviolation - Plot max. constraint violation at each iterationoptimplotfirstorderopt - Plot first-order optimality at each iterationoptimplotresnorm - Plot value of the norm of residuals at each iterationoptimplotstepsize - Plot step size at each iteration到此matlab中优化命令罗列完毕。,THANK YOU FOR YOUR PRESENCE,