Matlab优化工具箱简介.ppt

第五讲Matlab优化工具箱简介 optimizationtoobox 1 线性优化2 非线性优化3 极小化极大 Minmax 问题4 曲线拟合与插值 线性规划问题是目标函数和约束条件均为线性函数的问题 MATLAB7 0解决的线性规划问题的标准形式为minsub to 其中f x b beq lb ub为向量 A Aeq为矩阵 其它形式的线性规划问题都可经过适当变换化为此标准形式 在MATLAB5 x以上版中 线性规划问题LinearProgramming已用函数linprog取代了MATLAB5 x版中的lp函数 当然 由于版本的向下兼容性 一般说来 低版本中的函数在7 0版中仍可使用 5 1线性优化 函数linprog格式x linprog f A b 求minf xsub to线性规划的最优解 x linprog f A b Aeq beq 等式约束 若没有不等式约束 则A b x linprog f A b Aeq beq lb ub 指定x的范围 若没有等式约束 则Aeq beq x linprog f A b Aeq beq lb ub x0 设置初值x0 x linprog f A b Aeq beq lb ub x0 options options为指定的优化参数 x fval linprog 返回目标函数最优值 即fval f x x lambda exitflag linprog lambda为解x的Lagrange乘子 x lambda fval exitflag linprog exitflag为终止迭代的错误条件 x fval lambda exitflag output linprog output为关于优化的一些信息 说明 若exitflag 0表示函数收敛于解x exitflag 0表示超过函数估值或迭代的最大次数 exitflag 0表示函数不收敛于解x 若lambda lower表示下界lb lambda upper表示上界ub lambda ineqlin表示不等式约束 lambda eqlin表示等式约束 lambda中的非0元素表示对应的约束是有效约束 output iterations表示迭代次数 output algorithm表示使用的运算规则 output cgiterations表示PCG迭代次数 例5 1求下面的优化问题 min sub to 解 f 5 4 6 写成行向量亦可 A 1 11 324 320 b 20 42 30 lb zeros 3 1 x fval exitflag output lambda linprog f A b lb 结果为 x 最优解0 000015 00003 0000fval 最优值 78 0000exitflag 收敛1 output iterations 6 迭代次数cgiterations 0algorithm lipsol 所使用规则lambda ineqlin 3x1double eqlin 0 x1double upper 3x1double lower 3x1double lambda ineqlinans 0 00001 50000 5000 lambda lower ans 1 00000 00000 0000表明 不等约束条件2和3以及第1个下界是有效的 请写出下面线性规划的Matlab程序 c 0 4 0 28 0 32 0 72 0 64 0 6 A 0 010 010 010 030 030 03 0 02000 0500 00 02000 050 000 03000 08 b 850 700 100 900 Aeq beq vlb 0 0 0 0 0 0 vub x fval linprog c A b Aeq beq vlb vub MATLAB求解优化问题的主要函数 优化函数的输入变量 优化函数的输出变量 5 2非线性优化 5 2 1有约束的一元函数的最小值单变量函数求最小值的标准形式为sub to函数fminbnd格式x fminbnd fun x1 x2 返回自变量x在区间上函数fun取最小值时x值 fun为目标函数的表达式字符串或MATLAB自定义函数的函数柄 x fminbnd fun x1 x2 options x fval fminbnd x fval exitflag fminbnd x fval exitflag output fminbnd 例5 2计算下面函数在区间 0 1 内的最小值 解 x fval exitflag output fminbnd x 3 cos x x log x exp x 0 1 x 0 5223fval 0 3974exitflag 1 output iterations 9funcCount 9algorithm goldensectionsearch parabolicinterpolation 例5 3在 0 5 上求下面函数的最小值解 先自定义函数 在MATLAB编辑器中建立M文件为 functionf myfun x f x 3 2 1 保存为myfun m 然后在命令窗口键入命令 x fminbnd myfun 0 5 则结果显示为 x 3 5 2 2无约束多元函数最小值多元函数最小值的标准形式为其中 x为向量 命令利用函数fminsearch求无约束多元函数最小值 函数fminsearch格式x fminsearch fun x0 x0为初始点 fun为目标函数的表达式字符串或MATLAB自定义函数的函数柄 x fminsearch fun x0 options options查optimset x fval fminsearch 最优点的函数值 x fval exitflag fminsearch exitflag与单变量情形一致 x fval exitflag output fminsearch output与单变量情形一致 例5 4求的最小值点 解 X fminsearch 2 x 1 3 4 x 1 x 2 3 10 x 1 x 2 x 2 2 0 0 结果为X 1 00160 8335或在MATLAB编辑器中建立函数文件 functionf myfun x f 2 x 1 3 4 x 1 x 2 3 10 x 1 x 2 x 2 2 保存为myfun m 在命令窗口键入 X fminsearch myfun 0 0 或 X fminsearch myfun 0 0 结果为 X 1 00160 8335 5 2 3有约束的多元函数最小值非线性有约束的多元函数的标准形式为 sub to其中 x b beq lb ub是向量 A Aeq为矩阵 C x Ceq x 是返回向量的函数 f x 为目标函数 f x C x Ceq x 可以是非线性函数 在MATLAB5 x中 它的求解由函数constr实现 函数fmincon格式x fmincon fun x0 A b x fmincon fun x0 A b Aeq beq x fmincon fun x0 A b Aeq beq lb ub x fmincon fun x0 A b Aeq beq lb ub nonlcon x fmincon fun x0 A b Aeq beq lb ub nonlcon options x fval fmincon x fval exitflag fmincon x fval exitflag output fmincon x fval exitflag output lambda fmincon x fval exitflag output lambda grad fmincon x fval exitflag output lambda grad hessian fmincon 参数说明 fun为目标函数 它可用前面的方法定义 nonlcon的作用是通过接受的向量x来计算非线性不等约束和等式约束分别在x处的估计C和Ceq 通过指定函数柄来使用 如 x fmincon myfun x0 A b Aeq beq lb ub mycon 先建立非线性约束函数 并保存为mycon m function C Ceq mycon x C 计算x处的非线性不等约束的函数值 Ceq 计算x处的非线性等式约束的函数值 lambda是Lagrange乘子 它体现哪一个约束有效 output输出优化信息 grad表示目标函数在x处的梯度 hessian表示目标函数在x处的Hessian值 例5 5求下面问题在初始点 0 1 处的最优解 min sub to 解 约束条件的标准形式为 sub to先在MATLAB编辑器中建立非线性约束函数文件 function c ceq mycon x c x 1 1 2 x 2 ceq 无等式约束 然后 在命令窗口键入如下命令或建立M文件 fun x 1 2 x 2 2 x 1 x 2 2 x 1 5 x 2 目标函数 x0 01 A 23 线性不等式约束 b 6 Aeq 无线性等式约束 beq lb x没有下 上界 ub x fval exitflag output lambda grad hessian fmincon fun x0 A b Aeq beq lb ub mycon 则结果为x 34fval 13exitflag 1 解收敛 output iterations 2funcCount 9stepsize 1algorithm medium scale SQP Quasi Newton line search firstorderopt cgiterations lambda lower 2x1double x下界有效情况 通过lambda lower可查看 upper 2x1double x上界有效情况 为0表示约束无效 eqlin 0 x1double 线性等式约束有效情况 不为0表示约束有效 eqnonlin 0 x1double 非线性等式约束有效情况 ineqlin 2 5081e 008 线性不等式约束有效情况 neqnonlin 6 1938e 008 非线性不等式约束有效情况 grad 目标函数在最小值点的梯度 1 0e 006 0 1776hessian 目标函数在最小值点的Hessian值 1 0000 0 0000 0 00001 0000 5 2 4二次规划问题二次规划问题 quadraticprogramming 的标准形式为 sub to其中 H A Aeq为矩阵 f b beq lb ub x为向量其它形式的二次规划问题都可转化为标准形式 MATLAB5 x版中的qp函数已被6 0版中的函数quadprog取代 函数quadprog格式x quadprog H f A b 其中H f A b为标准形中的参数 x为目标函数的最小值 x quadprog H f A b Aeq beq Aeq beq满足等约束条件 x quadprog H f A b Aeq beq lb ub lb ub分别为解x的下界与上界 x quadprog H f A b Aeq beq lb ub x0 x0为设置的初值x quadprog H f A b Aeq beq lb ub x0 options options为指定的优化参数 x fval quadprog fval为目标函数最优值 x fval exitflag quadprog exitflag与线性规划中参数意义相同 x fval exitflag output quadprog output与线性规划中参数意义相同 x fval exitflag output lambda quadprog lambda与线性规划中参数意义相同 例5 6求二次规划的最优解maxf x1 x2 x1x2 3sub tox1 x2 2 0解 化成标准形式 sub tox1 x2 2 在Matlab中实现如下 H 0 1 1 0 f 0 0 Aeq 11 b 2 x fval exitflag output lambda quadprog H f Aeq b 结果为 x 1 00001 0000 fval 1 0000exitflag 4output iterations 1algorithm large scale projectivepreconditionedconjugategradients firstorderopt 0cgiterations 1message Optimizationterminated localminimumfound thesolutionissingular lambda eqlin 1 0000ineqlin lower upper 5 3极小化极大 Minmax 问题 sub to 其中 x b beq lb ub是向量 A Aeq为矩阵 C x Ceq x 和F x 是返回向量的函数 F x C x Ceq x 可以是非线性函数 函数fminimax格式x fminimax fun x0 x fminimax fun x0 A b x fminimax fun x0 A b Aeq beq x fminimax fun x0 A b Aeq beq lb ub x fminimax fun x0 A b Aeq beq lb ub nonlcon x fminimax fun x0 A b Aeq beq lb ub nonlcon options x fval maxfval fminimax x fval maxfval exitflag fminimax x fval maxfval exitflag output fminimax x fval maxfval exitflag output lambda fminimax 例5 7求下列函数最大值的最小化问题其中 解 先建立目标函数文件 并保存为myfun m functionf 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 然后 在命令窗口键入命令 x0 0 1 0 1 初始值 x fval fminimax myfun x0 结果为 x 4 00004 0000fval 0 0000 64 0000 2 0000 8 0000 0 0000 5 4曲线拟合与插值 在大量的应用领域中 人们经常面临用一个解析函数描述数据 通常是测量值 的任务 对这个问题有两种方法 插值 在插值法里 数据假定是正确的 要求以某种方法描述数据点之间所发生的情况 曲线拟合 曲线拟合或回归是人们设法找出某条光滑曲线 它最佳地拟合数据 但不必要经过任何数据点 标有 o 的是数据点 连接数据点的实线描绘了线性内插 虚线是数据的最佳拟合 曲线拟合的两个基本问题 1 最佳拟合意味着什么 2 应该用什么样的曲线 可用许多不同的方法定义最佳拟合 并存在无穷数目的曲线 当最佳拟合被解释为在数据点的最小误差平方和 且所用的曲线限定为多项式时 那么曲线拟合是相当简捷的 数学上 称为多项式的最小二乘曲线拟合 1Matlab曲线拟合 虚线和标志的数据点之间的垂直距离是在该点的误差 对各数据点距离求平方 并把平方距离全加起来 就是误差平方和 这条虚线是使误差平方和尽可能小的曲线 即是最佳拟合 最小二乘这个术语仅仅是使误差平方和最小的省略说法 Matlab曲线拟合和插值命令 x 0 1 2 3 4 5 6 7 8 91 y 4471 9783 286 167 087 347 669 569 489 3011 2 n 2 p polyfit x y n ezplot 9 8108 x x 20 1293 x 0 0317 二次多项式系数 既是p的输出 该命令为画出拟合曲线 xi linspace 0 1 100 x axisdataforplottingz polyval p xi 求出多项式在xi点处的取值 plot x y o x y xi z xlabel x ylabel y f x title SecondOrderCurveFitting 多项式阶次的选择是任意的 两点决定一直线或一阶多项式 三点决定一个平方或2阶多项式 按此进行 n 1数据点唯一地确定n阶多项式 于是 在上面的情况下 有11个数据点 我们可选一个高达10阶的多项式 然而 高阶多项式给出很差的数值特性 人们不应选择比所需的阶次高的多项式 此外 随着多项式阶次的提高 近似变得不够光滑 因为较高阶次多项式在变零前 可多次求导 原始数据标以 o 2阶曲线拟合是虚线 10阶拟合是实线 注意 在10阶拟合中 在左边和右边的极值处 数据点之间出现大的纹波 当企图进行高阶曲线拟合时 这种纹波现象经常发生 显然 越多就越好 的观念在这里不适用 2插值命令 插值定义为对数据点之间函数的估值方法 这些数据点是由某些集合给定 当人们不能很快地求出所需中间点的函数值时 插值是一个有价值的工具 例如 当数据点是某些实验测量的结果或是过长的计算过程时 就有这种情况 最简单插值的例子是MATLAB的作图 按缺省 MATLAB用直线连接所用的数据点以作图 这个线性插值猜测中间值落在数据点之间的直线上 当然 当数据点个数的增加和它们之间距离的减小时 线性插值就更精确 例如 x1 linspace 0 2 pi 60 x2 linspace 0 2 pi 6 plot x1 sin x1 x2 sin x2 xlabel x ylabel sin x title LinearInterpolation 一个在数据点之间用60个点 它比另一个只用6个点更光滑和更精确 如曲线拟合一样 插值要作决策 根据所作的假设 有多种插值 而且 可以在一维以上空间中进行插值 即如果有反映两个变量函数的插值 z f x y 那么就可在x之间和在y之间 找出z的中间值进行插值 MATLAB在一维函数interp1和在二维函数interp2中 提供了许多的插值选择 为了说明一维插值 考虑下列问题 12小时内 一小时测量一次室外温度 数据存储在两个MATLAB变量中 hours 1 12 indexforhourdatawasrecordedtemps 589152529313022252724 recordedtemperaturesplot hours temps hours temps viewtemperaturestitle Temperature xlabel Hour ylabel DegreesCelsius 正如上图看到的 MATLAB画出了数据点线性插值的直线 为了计算在任意给定时间的温度 人们可试着对可视的图作解释 另外一种方法 可用函数interp1 t interp1 hours temps 9 3 estimatetemperatureathour 9 3t 22 9000t interp1 hours temps 4 7 estimatetemperatureathour 4 7t 22t interp1 hours temps 3 26 57 111 7 findtempatmanypoints t 10 200030 000030 900024 9000 interp1的缺省用法是由interp1 x y xo 来描述 这里x是独立变量 横坐标 y是应变量 纵坐标 xo是进行插值的一个数值数组 另外 该缺省的使用假定为线性插值 若不采用直线连接数据点 我们可采用某些更光滑的曲线来拟合数据点 最常用的方法是用一个3阶多项式 即3次多项式 来对相继数据点之间的各段建模 每个3次多项式的头两个导数与该数据点相一致 这种类型的插值被称为3次样条或简称为样条 函数interp1也能执行3次样条插值 t interp1 hours temps 9 3 spline estimatetemperatureathour 9 3t 21 8577t interp1 hours temps 4 7 spline estimatetemperatureathour 4 7t 22 3143t interp1 hours temps 3 26 57 111 7 spline t 9 673430 042731 175525 3820 interp1二个强约束 1 人们不能要求有独立变量范围以外的结果 例如 interp1 hours temps 13 5 导致一个错误 因为hours在1到12之间变化 2 独立变量必须是单调的 即独立变量在值上必须总是增加的或总是减小的 二维插值是基于与一维插值同样的基本思想 然而 正如名字所隐含的 二维插值是对两变量的函数z f x y 进行插值 3 非线性数据 曲线 拟合非线性曲线拟合是已知输入向量xdata和输出向量ydata 并且知道输入与输出的函数关系为ydata F x xdata 但不知道系数向量x 今进行曲线拟合 求x使得下式成立 函数lsqcurvefit格式x lsqcurvefit fun x0 xdata ydata x lsqcurvefit fun x0 xdata ydata lb ub options x resnorm lsqcurvefit x resnorm residual exitflag output lambda lsqcurvefit x resnorm residual exitflag output lambda jacobian lsqcurvefit resnorm sum fun x xdata ydata 2 即在x处残差的平方和 residual fun x xdata ydata 即在x处的残差 exitflag为终止迭代的条件 output为输出的优化信息 lambda为解x处的Lagrange乘子 jacobian为解x处拟合函数fun的jacobian矩阵 例5 8求解如下最小二乘非线性拟合问题已知输入向量xdata和输出向量ydata 且长度都是n 拟合函数为即目标函数为其中 初始解向量为x0 0 3 0 4 0 1 解 先建立拟合函数文件 并保存为myfun mfunctionF myfun x xdata F x 1 xdata 2 x 2 sin xdata x 3 xdata 3 然后给出数据xdata和ydata xdata 3 67 79 34 18 62 81 37 910 05 4 ydata 16 5150 6263 124 7208 59 92 7163 9325 054 3 x0 10 10 10 初始估计值 x resnorm lsqcurvefit myfun x0 xdata ydata 结果为 Optimizationterminatedsuccessfully RelativefunctionvaluechangingbylessthanOPTIONS TolFunx 0 22690 33850 3021resnorm 6 2950 4非线性最小二乘非线性最小二乘 非线性数据拟合 的标准形式为 其中 L为常数设则目标函数可表达为其中 x为向量 F x 为函数向量 函数lsqnonlin格式x lsqnonlin fun x0 x0为初始解向量 fun为 i 1 2 m fun返回向量值F 而不是平方和值 平方和隐含在算法中 fun的定义与前面相同 x lsqnonlin fun x0 lb ub options options为指定优化参数 若x没有界 则lb ub x resnorm lsqnonlin resnorm sum fun x 2 即解x处目标函数值 x resnorm residual lsqnonlin residual fun x 即解x处fun的值 例5 9求下面非线性最小二乘问题初始解向量为x0 0 3 0 4 解 先建立函数文件 并保存为myfun m 由于lsqnonlin中的fun为向量形式而不是平方和形式 因此 myfun函数应由建立 k 1 2 10functionF myfun x k 1 10 F 2 2 k exp k x 1 exp k x 2 然后调用优化程序 x0 0 30 4 x resnorm lsqnonlin myfun x0 x 0 25780 2578resnorm 求目标函数值124 3622 实际问题求解 任务分配问题 某车间有甲 乙两台机床 可用于加工三种工件 假定这两台车床的可用台时数分别为800和900 三种工件的数量分别为400 600和500 且已知用三种不同车床加工单位数量不同工件所需的台时数和加工费用如下表 问怎样分配车床的加工任务 才能既满足加工工件的要求 又使加工费用最低 解设在甲车床上加工工件1 2 3的数量分别为x1 x2 x3 在乙车床上加工工件1 2 3的数量分别为x4 x5 x6 可建立以下线性规划模型 f 1391011128 A 0 41 110000000 51 21 3 b 800 900 Aeq 100100010010001001 beq 400600500 vlb zeros 6 1 vub x fval linprog f A b Aeq beq vlb vub 程序详解 见example5 m Theend