最优化方法概述

1、最优化方法概述,最优化技术是一门较新的学科分支。它是在上世纪五十年代初在电子计算机广泛应用的推动下才得到迅速发展，并成为一门直到目前仍然十分活跃的新兴学科。 最优化方法，顾名思义，就是从所有可能方案中选择最合理的一种以达到最优目标的方法。由于其宗旨是追求最优目标，因此它具有广泛的应用性。,最优化问题的一般形式,其中x= x1，x2，xnT,最优化问题的一般形式,目标函数,不等式约束,等式约束,x称为决策变量， 满足所有约束的变量称为可行解或可行点，可行点的集合称为可行域。 问题的求解是指在可行域中找一点x*，使得目标函数在该点取极小值，这样的点称为问题的最优点，也称为最小点，而相应的目标函数值

2、f(x*)称为最优值，(x*,f(x*)称为最优解，习惯上x*称为最优解。,最优化问题常用的概念,定义1：整体（全局）最优解：若 ，对于一切 ，恒有 则称 是最优化问题的整体最优解。 定义2：局部最优解：若 ，存在某邻域 ，使得对于一切 ，恒有 则称 是最优化问题的局部最优解。其中 严格最优解：当 ，有 则称 为问题的严格最优解。,f(X),局部最优解,整体最优解,最优化问题分类,最优化方法的基本方法迭代法,为什么求最优解要用迭代法？下面举例说明。 例 则最优解必定是方程 的解 ， 但x 无法显化。而由于 可知 的解存在,最优化方法的基本方法迭代法,下降：在每次迭代中，后继点处的函数值要有所减

3、少。,迭代：,最优化方法的基本方法迭代法,迭代算法的步骤：,线性规划是求线性目标函数在线性约束条件下的最大值或最小值的问题。,线性规划问题的数学模型是将实际问题转化为一组线性不等式或等式约束下求线性目标函数的最小(大)值问题, 它都可以化为如下标准(矩阵)形式：,线性规划,单纯形法,单纯形方法基本思路： 从可行域中某个基础可行解（一个顶点）开始（称为初始基础可行解）。 如可能，从可行域中求出具有更优目标函数值的另一个基础可行解（另一个顶点），以改进初始解。 继续寻找更优的基础可行解，进一步改进目标函数值。当某一个基础可行解不能再改善时，该解就是最优解。,数学模型 max S = 50 x1 +

4、 30 x2 s.t. 4x1 + 3x2 120 2x1 + x2 50 x1 , x2 0,数学模型的标准形: maxS = 50 x1 + 30 x2 s.t. 4x1 + 3x2+ x3 = 120 2x1+ x2 + x4= 50 x1,x2 , x3 , x4 0,基础可行解X(1)=(0, 0,120,50)T, 目标函数值 S= 0,X(1)=(0, 0,120,50)T 相当于O(0,0),X1 进基变量， X4出基变量， 主元（2）,第二行除以2,第一行加上第二行的负4倍,新的基础可行解X(2)=(25, 0,20,0)t, 目标函数值 S = 1250 计算检验数,X(2

5、)=(25, 0,20,0)t 相当于Q1(25,0),X2 进基变量， X3出基变量， 主元（1）,第二行加上第一行的负1/2倍,基础可行解X(3)=(15,20,0,0)T, 也是最优解,其最优值 S = 1350,X(3)=(15,20,0,0)T 相当于Q2(15,20),X(1)=(0, 0,120,50)T 相当于O(0,0),X(1)=(0, 0,120,50)T 相当于O(0,0),X(2)=(25, 0,20,0)T 相当于Q1(25,0),X(2)=(25, 0,20,0)T 相当于Q1(25,0),X(3)=(15,20,0,0)T 相当于Q2(15,20),用Matla

6、b求解线性规划,模型 用命令 x, fval= linprog(f,A,b,A1,b1,lb,ub),minf(x)=-5x1-4x2-6x3 s.t. x1-x2+x320 3x1+2x2+4x342 3x1+2x230 0x1, 0x2,0x3,用Matlab求解线性规划,用Matlab求解线性规划,解问题 把问题极小化并将约束标准化,用Matlab求解线性规划,键入c=-2,-3,5;a=-2,5,-1; b=-10;a1=1,1,1;b1=7;LB=0,0,0; x,y=linprog(c,a,b,a1,b1,LB) 得当X=(6.4286,0.5714,0.0000)时, z=-14

7、.5714最大.,用Matlab求解线性规划,无约束非线性规划,min f (x)| x En , 这里x =(x1 , x2 , , xn)T.,一元函数无约束优化问题,多元函数无约束优化问题,精确线搜索 试探法: 黄金分割法、Fibonacci法、二分法 函数逼近法: Newton法、割线法、抛物线法、 三次插值法 非精确线搜索 Armijo步长规则、Goldstein步长规则、 Wolfe步长规则,一元函数无约束优化问题(一维搜索),分割法原理：设函数 f (x) 在闭区间 a, b 上是下单峰函数, 即在 (a, b) 内 f (x) 由唯一的极小点x*, 在x* 的左边 f (x)

8、严格单调下降, 在x* 的右边f (x)严格单调上升. 那么对于(a, b)内任意两点x1x2, 如果 f (x1) f (x2 ), 则x*a, x2；否则x*x1, b.,如右图,我们先介绍求解一元函数 y = f (x)极小值的数值迭代算法：一维搜索算法中的黄金分割法(0.618法).,一元函数无约束优化问题,取 x2 = a + 0.618 (b - a), f2 = f (x2), 转向. 取 x1 = a + 0.382 (b - a), f1 = f (x1), 转向. 若 | b a | , 则取x* = (a + b )/2, 停. 否则转向. 若f1 f2 , 则取b =

9、x2 , x2 = x1, f2 = f1 , 转向; 若f1= f2 , 则取a = x1, b = x2, 转向; 若f1f2 , 则取a = x1, x1= x2, f1 = f2 , 转向. 取 x2 = a + 0.618 (b - a), f2= f (x2), 转向.,给定下单峰区间 a, b 及控制误差0, 黄金分割法(0.618法)的迭代步骤：,函数逼近法：牛顿法,基本思想：在极小点附近用二阶Taylor多项式近似。,算法步骤:,缺点：初始点选择十分重要。如果初始点靠近极小点，则 可能很快收敛；如果初始点远离极小点，迭代产生的点 列可能不收敛于极小点。,例：,解：,用Matl

10、ab解无约束优化问题,其中（3）、（4）、（5）的等式右边可选用（1）或（2）的等式右边。 函数fminbnd的算法基于黄金分割法，它要求目标函数必须是连续函数，并可能只给出局部最优解。,常用格式如下： （1）x= fminbnd (fun,x1,x2) （2）x= fminbnd (fun,x1,x2 ，options) （3）x，fval= fminbnd（.） （4）x，fval，exitflag= fminbnd（.） （5）x，fval，exitflag，output= fminbnd（.）,主程序为wliti1.m: f=2*exp(-x).*sin(x); fplot(f,0,8

11、); %作图语句 xmin,ymin=fminbnd (f, 0,8) f1=-2*exp(-x).*sin(x); xmax,ymax=fminbnd (f1, 0,8),例2 对边长为3米的正方形铁板，在四个角剪去相等的正方形以制成方形无盖水槽，问如何剪法使水槽的容积最大？,解,先编写M文件fun0.m如下: function f=fun0(x) f=-(3-2*x).2*x;,主程序为wliti2.m: x,fval=fminbnd(fun0,0,1.5); xmax=x fmax=-fval,运算结果为: xmax = 0.5000,fmax =2.0000.即剪掉的正方形的边长为0.

12、5米时水槽的容积最大,最大容积为2立方米.,问题：,在点,处，,沿什么方向,下降最快?,结论：,负梯度方向使,下降最快，,亦即最速,下降方向,最速下降法,最速下降法算法,Step1:,给出,Step2:,计算,如果,停,Step3:,计算下降方向,Step4:,计算步长因子,Step5:,令,转步.,分析：,设,是正定二次函数,由精确的线搜索确定的,特别当：,例1：,用最速下降法求解：,解：,分析：,(1),因此：,最速下降法是整体收敛的，,且是线性收敛的,(2),两个相邻的搜索方向是正交的,最速下降法优点,(1),程序设计简单，计算量小，存储量小，,对初始点没有特别要求,(2),有着很好的整

13、体收敛性，即使对一般的,目标函数，它也整体收敛,最速下降法缺点,(1),最速下降法是线性收敛的，并且有时是,很慢的线性收敛,原因：,仅反映,在,处,的局部性质,相继两次迭代中搜索,方向是正交的,小结,(1),最速下降法是基本算法之一，而非有效,的实用算法,最速下降法的本质是用线性函数来近似,目标函数，,要想得到快速算法，需要考,虑对目标函数的高阶逼近,基本思想,利用目标函数,在点,处的二阶Taylor,展开式去近似目标函数，,用二次函数的极小点,去逼近目标函数的极小点,牛顿法,算法构造,问题：,设,二阶连续可微，,海赛阵,正定,如何从,因为,正定，,则,有唯一极小点，,用这个,极小点作为,所以

14、要求：,即：,因此：,这就是牛顿法迭代公式,注：,这里,牛顿法算法,Step1:,给出,Step2:,计算,如果,停,Step3:,否则计算,Step4:,令,转步.,并且求解方程,得出,例1：,用牛顿法求解：,解：,牛顿法优点,(1),(2),对正定二次函数，迭代一次就可以得到,极小点,如果,正定且初始点选取合适，,算法,二阶收敛,牛顿法缺点,(1),(2),对多数问题算法不是整体收敛的,每次都需要计算,计算量大,(3),每次都需要解,方程组有时奇异或病态的，,无法确定,或,不是下降方向,(4),收敛到鞍点或极大点的可能性并不小,阻尼牛顿法算法,Step1:,给出,Step2:,计算,如果,

15、停,Step3:,否则计算,Step4:,沿,并且求解方程,得出,进行线搜索，,得出,Step5:,令,转Step2.,可以克服牛顿法的缺点（1）,共轭方向法 与共轭梯度法,算法特点,（1）有效的算法，克服了最速下降法的慢,收敛性，又避免了牛顿法的计算量大和局部收,性的缺点,（2）算法简单，易于编程，需存储空间小等,优点，是求解大规模问题的主要方法,共轭方向及其性质,定义1:,设,是,中任一组,非零向量，,如果：,则称,是关于,共轭的,注：,若,则是正交的，因此共轭是,正交的推广,共轭方向法算法,Step1:,给出,计算,和初始下降方向,Step2:,如果,停止迭代,Step3:,计算,使得,

16、Step4:,采用某种共轭方向法计算,使得：,Step5:,令,转Step2.,共轭梯度法,记：,左乘,并使,得：,（Hestenes-Stiefel公式）,取：,共轭梯度法基本性质,定理3:,对于正定二次函数，,采用精确线搜索,的共轭梯度法在,步后终止，,且对,成立下列关系式：,（共轭性）,（正交性）,（下降条件）,系数的其他形式,（）FR公式,（1964）,（2）PRP公式,（1969）,FR共轭梯度法算法,Step1:,给出,Step2:,计算,如果,停,Step3:,Step4:,由精确线搜索求,Step5:,转Step2.,例4：,用FR共轭梯度法求解：,解：,化成,形式,(1),(

17、2),例5：,用FR共轭梯度法求解：,解：,化成,形式,(1),(2),所用函数,x向量； f(x)返回值为标量的函数； fminunc寻找无约束多变量函数f(x)最小值.,用Matlab求解多元函数无约束优化问题,x = fminunc(fun,x0) 初始点为x0，寻找fun所定义的函数的局部最小。 x = fminunc(fun,x0,options) 初始点为x0，寻找fun所定义的函数的局部最小，优化选项包含在结构options中。 x,fval = fminunc(.) x为最小点，fval为最优值。 x,fval,exitflag = fminunc(.) exitflag标识算

18、法终止的原因。,x,fval,exitflag,output = fminunc(.) 结构output中包含优化信息。 x,fval,exitflag,output,grad = fminunc(.) 在grad中包含函数fun在x处的梯度值。 x,fval,exitflag,output,grad,hessian = fminunc(.) hessian中返回函数fun在x处的Hessian矩阵值。,问题,Matlab工具箱求解过程 步骤1: 写一个目标函数的M文件objfun.m function f = objfun(x) f = exp(x(1)*(4*x(1)2+2*x(2)2+4

19、*x(1)*x(2)+2*x(2)+1);,步骤2: 调用一个无约束优化程序 x0 = -1,1; % 起始点 options = optimset(LargeScale,off); x,fval,exitflag,output = fminunc(objfun,x0,options) 注：1.可把步骤二的内容建立一个M文件Unconstr.m,在命令行运行Unconstr即得到以下运行结果。 2.以M文件函数给出时文件名前加,运行结果 Matlab命令窗口显示如下运行结果。 Optimization terminated: relative infinity-norm of gradient

20、 less than options.TolFun. x = 0.5000 -1.0000 fval = 1.0983e-015 exitflag = 1 output = iterations: 8 funcCount: 66 stepsize: 1 firstorderopt: 7.3704e-008 algorithm: medium-scale: Quasi-Newton line search message: 1x85 char,结果分析 66次函数计算后，找到最终答案。 x = 0.5000 -1.0000 在终点x 处函数值为fval： fval = 1.0983e-015 e

21、xitflag 是算法是否收敛的标识. exitflag 0 意味着找到局部极小。 exitflag = 1 输出结构给出更多关于优化的信息。对于fminunc，包括iterations中的迭代数，funcCount中的函数计算次数，stepsize 中的最终步长，firstorderopt中的一阶最优性的测量。Algorithm中的算法的类型。,当多于一种局部最小存在，初始点 x1，x2的设置影响函 数评价的次数和最终点的值。在前面的例子中， x0被初始化为 -1，1。 变量options可被传递到 fminunc改变优化算法的特性， x = fminunc(objfun,x0,option

22、s)； option是一种包含终止误差和算法选择的结构。Options可以使用optimset函数创建。 options = optimset(LargeScale,off)； 这个例子中，已关闭了large-scale算法缺省选项，而使用medium-scale算法。其他选项包括控制优化迭代过程命令行显示的数量，终止准则的公差，使用使用者提供的梯度还是雅可比矩阵，迭代或函数计算的最大数。,非线性不等式约束,所用函数： 式中： x，b，beq，lb，ub 向量； A，Aeq矩阵； c(x)，ceq(x)返回向量的函数； f(x)返回标量的函数； f(x)，c(x)，ceq(x)可以是非线性函数

23、。,x = fmincon(fun,x0,A,b) x0为起始点，寻找函数fun的最小值，约束条件为线性不等式约束A*xb。x0可以是标量、向量或矩阵。 x = fmincon(fun,x0,A,b,Aeq,beq) 寻找函数fun的最小值，约束条件为线性等式约束Aeq*x = beq以及线性不等式约束A*x b。如果无线性不等式约束设置A=、b=。 x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub) 定义了变量x的上下界，使最优解x为lb x ub。如果等式约束不存在，设置Aeq=、beq= 如果x(i)无下界，设置lb(i)=-Inf；如果x(i) 无上界，设置ub

24、(i)=Inf。,x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon) nonlcon中定义了非线性不等式c(x)或等式ceq(x)约束。优化结果使c(x) 0、ceq(x) = 0。如果变量x无上下界，设置lb=、ub=。 x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options) 优化选项在结构options中定义。使用optimset 函数定义options. 如果没有非线性不等式和等式约束，设置nonlcon = 。 x,fval = fmincon(.) 返回最优解（x，fval）。 x,fval,

25、exitflag = fmincon(.) 返回描述fmincon 结束的原因的标识符exitflag。,x,fval,exitflag,output = fmincon(.) 结构output中包含优化信息。 x,fval,exitflag,output,lambda = fmincon(.) 返回最优解x的拉格朗日乘子，存储在结构lambda中。 x,fval,exitflag,output,lambda,grad = fmincon(.) 在grad中包含函数fun在x处的梯度值。 x,fval,exitflag,output,lambda,grad,hessian = fmincon(

26、.) 在hessian中包含函数fun在x处的Hessian矩阵值。,问题,约束条件：,Matlab工具箱求解过程 不能在命令行把非线性的约束传给fmincon。但可创建一个M-file，confun.m，这个文件返回当前x的两个约束条件的变量c。接着调用带约束的优化函数fmincon。由于fmincon要求约束条件被写成,的形式，因此必须重写以下形式的约束。,步骤1：写一个约束条件的confun.m 文件。 function c, ceq = confun(x) c = 1.5 + x(1)*x(2) - x(1) - x(2); -x(1)*x(2) - 10; %非线性不等式约束 ceq

27、 = ; % 非线性等式约束,步骤2：写一个目标函数的M文件objfun.m function f = objfun(x) f = exp(x(1)*(4*x(1)2+2*x(2)2+4*x(1)*x(2)+2*x(2)+1); 步骤3：调用带约束的优化程序。 x0 = -1,1; % 起始点 options = optimset(LargeScale,off,Display,iter); x,fval,exitflag,output= fmincon(objfun,x0,confun,options),运行结果 x = -9.5474 1.0474 fval = 0.0236 exitfla

28、g = 1 output = iterations: 8 funcCount: 36 stepsize: 1 algorithm: medium-scale: SQP, Quasi-Newton, line-search firstorderopt: 8.5125e-007 cgiterations: message: 1x144 char,运行结果分析 36 次函数调用后，最优解为 x = -9.5474 1.0474 fval = 0.0236 也可计算最优解处的约束。 c,ceq = confun(x) 返回 c= 1.0e-14 * 0.1110 -0.1776 ceq = 可注意到约

29、束值小于等于0，即x满足。,带边界的约束,问题,约束条件：,Matlab工具箱求解过程 通过对约束优化函数定义简单的边界约束，变量x可以被限制在某一范围。对于函数fmincon，命令 x = fmincon(objfun,x0,lb,ub,confun,options); 限制 x 为 lb x ub。,步骤1：写一个约束条件的confun.m 文件。 function c, ceq = confun(x) c = 1.5 + x(1)*x(2) - x(1) - x(2); -x(1)*x(2) - 10; %非线性不等式约束 ceq = ; % 非线性等式约束 步骤2：写一个目标函数的M文

30、件objfun.m function f = objfun(x) f = exp(x(1)*(4*x(1)2+2*x(2)2+4*x(1)*x(2)+2*x(2)+1);,步骤3：调用带约束的优化程序 为使方程中x大于0，使用如下命令。 x0 = -1,1; % 估计问题的起始点 lb = 0,0; % 设置下界 ub = ; % 无上界 options = optimset(LargeScale,off,Display,iter); x,fval,exitflag,output= fmincon(objfun,x0,lb,ub,confun,options) c, ceq = confun(

31、x) 注意：为了向fmincon传递下界作为第七个变量，必须定义第三个到第六个变量的值。本例中因为没有线性不等式或线性等式约束，定义这些变量为空。,运行结果,x = 0 1.5000 fval = 8.5000 exitflag = 1 output = iterations: 3 funcCount: 15 stepsize: 1 algorithm: medium-scale: SQP, Quasi-Newton, line-search firstorderopt: 1.6307e-012 cgiterations: message: 1x144 char c = 0 -10 ceq = ,运行结果分析 当 lb或 ub 比 x包含的元素少，仅仅x中对应的元素被约束。如果仅仅变量的一部分被约束，在 lb中使用-inf作为无约束的下边界向量，在 ub中使用 inf 作为无约束的上变量。例如： lb = -inf 0； ub = 10 inf； 边界x110、x20，