Matlab的优化工具箱的几个应用函数及例子

Matlab的优化工具箱的几个应用函数及例子

利用Matlab的优化工具箱，可以求解线性规划、非线性规划和多目标规划问题。具体而言，包括线性、非线性最小化，最大最小化，二次规划，半无限问题，线性、非线性方程（组）的求解，线性、非线性的最小二乘问题。另外，该工具箱还提供了线性、非线性最小化，方程求解，曲线拟合，二次规划等问题中大型课题的求解方法，为优化方法在工程中的实际应用提供了更方便快捷的途径。9.1.1 优化工具箱中的函数

优化工具箱中的函数包括下面几类：

1．最小化函数

表9-1 最小化函数表

函 数

fgoalattain fminbnd fmincon fminimax

fseminf linprog quadprog

描 述

多目标达到问题 有边界的标量非线性最小化 有约束的非线性最小化 最大最小化 半无限问题 线性课题 二次课题

fminsearch, fminunc 无约束非线性最小化

2．方程求解函数

表9-2 方程求解函数表

函 数

\\ fsolve fzero

线性方程求解 非线性方程求解 标量非线性方程求解

描 述

3．最小二乘（曲线拟合）函数

表9-3 最小二乘函数表

函 描

数

\\ lsqlin lsqcurvefit lsqnonlin lsqnonneg

线性最小二乘 有约束线性最小二乘 非线性曲线拟合 非线性最小二乘 非负线性最小二乘

述

4．实用函数

表9-4 实用函数表

函 数

optimset 设置参数 optimget

5．大型方法的演示函数

描 述

表9-5 大型方法的演示函数表

函 描 数 述

circustent 马戏团帐篷问题—二次课题 molecule 用无约束非线性最小化进行分子组成求解 optdeblur 用有边界线性最小二乘法进行图形处理

6．中型方法的演示函数

表9-6 中型方法的演示函数表

函 数

bandemo dfildemo goaldemo optdemo tutdemo

香蕉函数的最小化 过滤器设计的有限精度 目标达到举例 演示过程菜单 教程演示

描 述

9.1.3 参数设置

利用optimset函数，可以创建和编辑参数结构；利用optimget函数，可以获得options优化参数。 ● optimget函数

功能：获得options优化参数。 语法：

val = optimget(options,'param')

val = optimget(options,'param',default) 描述：

val = optimget(options,'param') 返回优化参数options中指定的参数的值。只需要用参数开头的字母来定义参数就行了。 val = optimget(options,'param',default) 若options结构参数中没有定义指定参数，则返回缺省值。注意，这种形式的函数主要用于其它优化函数。 举例： 1．

下面的命令行将显示优化参数options返回到my_options结构中：

val = optimget(my_options,'Display') 2．

下面的命令行返回显示优化参数options到my_options结构中（就象前面的例子一样），但如果显示参数没有定义，则返回值'final':

optnew = optimget(my_options,'Display','final'); 参见： optimset

● optimset函数

功能：创建或编辑优化选项参数结构。 语法：

options = optimset('param1',value1,'param2',value2,...) optimset

options = optimset

options = optimset(optimfun)

options = optimset(oldopts,'param1',value1,...) options = optimset(oldopts,newopts) 描述：

options = optimset('param1',value1,'param2',value2,...) 创建一个称为options的优化选项参数，其中指定的参数具有指定值。所有未指定的参数都设置为空矩阵[]（将参数设置为[]表示当options传递给优化函数时给参数赋缺省值）。赋值时只要输入参数前面的字母就行了。

optimset函数没有输入输出变量时，将显示一张完整的带有有效值的参数列表。 options = optimset (with no input arguments) 创建一个选项结构options，

其中所有的元素被设置为[]。

options = optimset(optimfun) 创建一个含有所有参数名和与优化函数optimfun相关的缺省值的选项结构options。

options = optimset(oldopts,'param1',value1,...) 创建一个oldopts的拷贝，用指定的数值修改参数。 options = optimset(oldopts,newopts) 将已经存在的选项结构oldopts与新的选项结构newopts进行合并。newopts参数中的所有元素将覆盖oldopts参数中的所有对应元素。 举例：

1．下面的语句创建一个称为options的优化选项结构，其中显示参数设为'iter'，TolFun参数设置为1e-8:

options = optimset('Display','iter','TolFun',1e-8)

2．下面的语句创建一个称为options的优化结构的拷贝，改变TolX参数的值，将新值保存到optnew参数中:

optnew = optimset(options,'TolX',1e-4);

3．下面的语句返回options优化结构，其中包含所有的参数名和与fminbnd函数相关的缺省值：

options = optimset('fminbnd') 4．若只希望看到fminbnd函数的缺省值，只需要简单地键入下面的语句就行了：

optimset fminbnd

或者输入下面的命令，其效果与上面的相同：

optimset('fminbnd') 参见： optimget

9.1.4 模型输入时需要注意的问题 使用优化工具箱时，由于优化函数要求目标函数和约束条件满足一定的格式，所以需要用户在进行模型输入时注意以下几个问题： 1.目标函数最小化

优化函数fminbnd、fminsearch、fminunc、fmincon、fgoalattain、fminmax和lsqnonlin都要求目标函数最小化，如果优化问题要求目标函数最大化，可以通过使该目标函数的负值最小化即-f(x)最小化来实现。近似地，对于quadprog

函数提供-H和-f，对于linprog函数提供-f。 2.约束非正

优化工具箱要求非线性不等式约束的形式为Ci(x)≤0，通过对不等式取负可以达到使大于零的约束形式变为小于零的不等式约束形式的目的，如Ci(x)≥0形式的约束等价于- Ci(x)≤0；Ci(x)≥b形式的约束等价于- Ci(x)+b≤0。 3.避免使用全局变量

9.1.5 @（函数句柄）函数

MATLAB6.0中可以用@函数进行函数调用。@函数返回指定MATLAB函数的句柄，其调用格式为：

handle = @function

利用@函数进行函数调用有下面几点好处：

●

用句柄将一个函数传递给另一个函数；

●

减少定义函数的文件个数；

●

改进重复操作；

●

保证函数计算的可靠性。

下面的例子为humps函数创建一个函数句柄，并将它指定为fhandle变量。

fhandle = @humps;

同样传递句柄给另一个函数，也将传递所有变量。本例将刚刚创建的函数句柄传递给fminbnd函数，然后在区间[0.3,1]上进行最小化。 x = fminbnd (@humps, 0.3, 1) x =

0.6370

9.2 最小化问题

9.2.1 单变量最小化 9.2.1.1 基本数学原理

本节讨论只有一个变量时的最小化问题，即一维搜索问题。该问题在某些情况下可以直接用于求解实际问题，但大多数情况下它是作为多变量最优化方法的基础在应用，因为进行多变量最优化要用到一维搜索法。该问题的数学模型为：

x = fminbnd(inline('sin(x*x)'),x0)

同样，fun参数可以是一个包含函数名的字符串。对应的函数可以是M文件、内部函数或MEX文件。若fun='myfun'，则M文件函数myfun.m必须有下面的形式：

function f = myfun(x) f = ...

%计算x处的函数值。

若fun函数的梯度可以算得，且options.GradObj设为'on'（用下式设定）,

options = optimset('GradObj','on')

则fun函数必须返回解x处的梯度向量g到第二个输出变量中去。注意，当被调用的fun函数只需要一个输出变量时（如算法只需要目标函数的值而不需要其梯度值时），可以通过核对nargout的值来避免计算梯度值。 function [f,g] = myfun(x) f = ...

%计算x处得函数值。 if nargout > 1

%调用fun函数并要求有两个输出变量。 g = ...

%计算x处的梯度值。 end

若Hessian矩阵也可以求得，并且options.Hessian设为'on',即,

options = optimset('Hessian','on')

则fun函数必须返回解x处的Hessian对称矩阵H到第三个输出变量中去。注意，当被调用的fun函数只需要一个或两个输出变量时（如算法只需要目标函数的值f和梯度值g而不需要Hessian矩阵H时），可以通过核对nargout的值来避免计算Hessian矩阵 function [f,g,H] = myfun(x) f = ...

% 计算x处得函数值。 if nargout > 1

% 调用fun函数并要求有两个输出变量。 g = ...

% 计算x处的梯度值。

if nargout > 2 H = ...

% 计算x处的Hessian矩阵。 end

优化参数选项。可以通过optimset函数设置或改变这些参数。其中有的参数适用于所有的优化算法，有的则只适用于大型优化问题，另外一些则只适用于中型问题。

首先描述适用于大型问题的选项。这仅仅是一个参考，因为使用大型问题算法有一些条件。对于fminunc函数来说，必须提供梯度信息。 ??

LargeScale – 当设为'on'时使用大型算法，若设为'off'则使用中型问题的算法。 适用于大型和中型算法的参数： ??

Diagnostics – 打印最小化函数的诊断信息。 ??

Display – 显示水平。选择'off'，不显示输出；选择'iter'，显示每一步迭代过程的输出；选择'final'，显示最终结果。打印最小化函数的诊断信息。 ??

GradObj – 用户定义的目标函数的梯度。对于大型问题此参数是必选的，对于中型问题则是可选项。 ??

MaxFunEvals – 函数评价的最大次数。 ??

MaxIter – 最大允许迭代次数。 ??

TolFun – 函数值的终止容限。

options

??

TolX – x处的终止容限。 只用于大型算法的参数： ??

Hessian – 用户定义的目标函数的Hessian矩阵。 ??

HessPattern – 用于有限差分的Hessian矩阵的稀疏形式。若不方便求fun函数的稀疏Hessian矩阵H，可以通过用梯度的有限差分获得的H的稀疏结构（如非零值的位置等）来得到近似的Hessian矩阵H。若连矩阵的稀疏结构都不知道，则可以将HessPattern设为密集矩阵，在每一次迭代过程中，都将进行密集矩阵的有限差分近似（这是缺省设置）。这将非常麻烦，所以花一些力气得到Hessian矩阵的稀疏结构还是值得的。 ??

MaxPCGIter – PCG迭代的最大次数。 ??

PrecondBandWidth – PCG前处理的上带宽，缺省时为零。对于有些问题，增加带宽可以减少迭代次数。 ??

TolPCG – PCG迭代的终止容限。 ??

TypicalX – 典型x值。 只用于中型算法的参数：

??

DerivativeCheck – 对用户提供的导数和有限差分求出的导数进行对比。 ??

DiffMaxChange – 变量有限差分梯度的最大变化。 ??

DiffMinChange - 变量有限差分梯度的最小变化。 ??

LineSearchType – 一维搜索算法的选择。 描述退出条件: ??>0

表示目标函数收敛于解x处。

exitflag

??0

表示已经达到函数评价或迭代的最大次数。 ??<0

表示目标函数不收敛。 该参数包含下列优化信息： ??

output.iterations – 迭代次数。 ??

output.algorithm – 所采用的算法。 ??

output output.funcCount – 函数评价次数。

??

output.cgiterations – PCG迭代次数（只适用于大型规划问题）。 ??

output.stepsize – 最终步长的大小（只用于中型问题）。 ??

output.firstorderopt – 一阶优化的度量：解x处梯度的范数。

注意：

1．对于求解平方和的问题，fminunc函数不是最好的选择，用lsqnonlin函数效果更佳。 2．使用大型方法时，必须通过将options.GradObj设置为'on'来提供梯度信息，否则将给出警告信息。 算法：

大型优化算法

若用户在fun函数中提供梯度信息，则缺省时函数将选择大型优化算法，该算法是基于内部映射牛顿法的子空间置信域法，理论描述可参见文献[8],[9]。计算中的每一次迭代涉及到用PCG法求解大型线性系统得到的近似解。 中型优化算法

此时fminunc函数的参数options.LargeScale设置为'off'。该算法采用的是基于二次和三次混合插值一维搜索法的BFGS拟牛顿法。该法通过BFGS公式来更新Hessian矩阵。通过将HessUpdate参数设置为'dfp'，可以用DFP公式来求得Hessian矩阵逆的近似。通过将HessUpdate参数设置为'steepdesc'，可以用最速下降法来更新Hessian矩阵。但一般不建议使用最速下降法。

缺省时的一维搜索算法，当options.LineSearchType设置为'quadcubic'时,将采用二次和三次混合插值法。将options.LineSearchType设置为'cubicpoly'时，将采用三次插值法。第二种方法需要的目标函数计算次数更少，但梯度的计算次数更多。这样，如果提供了梯度信息，或者能较容易地算得，则三次插值法是更佳的选择。 局限性： 1．

目标函数必须是连续的。fminunc函数有时会给出局部最优解。 2．

fminunc函数只对实数进行优化，即x必须为实数，而且f(x)必须返回实数。当x为复数时，必须将它分解为实部和虚部。 3．

在使用大型算法时，用户必须在fun函数中提供梯度（options参数中GradObj属性必须设置为'on'）。 4． 目前，若在fun函数中提供了解析梯度，则options参数DerivativeCheck不能用于大型算法以比较解析梯度和有限差分梯度。通过将options参数的MaxIter 属性设置为0来用中型方法核对导数。然后重新用大型方法求解问题。 参见：

fminsearch, optimset, inline

fminsearch函数

功能：求解多变量无约束函数的最小值。 数学模型：

其中x为向量，f(x)为函数，返回标量。 语法：

x = fminsearch(fun,x0)

x = fminsearch(fun,x0,options)

x = fminsearch(fun,x0,options,P1,P2,...) [x,fval] = fminsearch(...)

[x,fval,exitflag] = fminsearch(...)

[x,fval,exitflag,output] = fminsearch(...) 描述：

fminsearch 求解多变量无约束函数的最小值。该函数常用于无约束非线性最优化问题。

x = fminsearch(fun,x0) 初值为x0，求fun函数的局部极小点x。x0可以是标量、向量或矩阵。

x = fminsearch(fun,x0,options)用options参数指定的优化参数进行最小化。 x = fminsearch(fun,x0,options,P1,P2,...) 将问题参数p1、p2等直接输给目标函数fun，将options参数设置为空矩阵，作为options参数的缺省值。 [x,fval] = fminsearch(...)将x处的目标函数值返回到fval参数中。 [x,fval,exitflag] = fminsearch(...)返回exitflag值，描述函数的退出条件。 [x,fval,exitflag,output] = fminsearch(...)返回包含优化信息的输出参数output。 变量：

各变量的意义同前。 算法：

fminsearch使用单纯形法进行计算。

对于求解二次以上的问题，fminsearch函数比fminunc函数有效。但是，当问题为高度非线性时，fminsearch函数更具稳健性。 局限性： 1．

应用fminsearch函数可能会得到局部最优解。 2．

fminsearch函数只对实数进行最小化，即x必须由实数组成，f(x)函数必须返回实数。如果x时复数，必须将它分为实数部和虚数部两部分。 注意:

fminsearch函数不适合求解平方和问题，用lsqnonlin函数更好一些。 参见：

fminbnd, fminunc, optimset, inline

9.2.4 二次规划问题 9.2.4.1 基本数学原理

如果某非线性规划的目标函数为自变量的二次函数，约束条件全是线性函数，就称这种规划为二次规划。其数学模型为：

其中，H, A,和Aeq为矩阵,f, b, beq, lb, ub,和x为向量。 9.2.4.2 相关函数介绍 quadprog函数

功能：求解二次规划问题。 语法：

x = quadprog(H,f,A,b)

x = quadprog(H,f,A,b,Aeq,beq)

x = quadprog(H,f,A,b,Aeq,beq,lb,ub) x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0)

x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options) [x,fval] = quadprog(...)

[x,fval,exitflag] = quadprog(...)

[x,fval,exitflag,output] = quadprog(...)

[x,fval,exitflag,output,lambda] = quadprog(...)

描述：

x = quadprog(H,f,A,b) 返回向量x，最小化函数 1/2*x'*H*x + f'*x ，其约束条件为 A*x <= b。

x = quadprog(H,f,A,b,Aeq,beq)仍然求解上面的问题，但添加了等式约束条件Aeq*x = beq。

x = quadprog(H,f,A,b,lb,ub)定义设计变量的下界lb和上界ub，使得lb <= x <= ub。

x = quadprog(H,f,A,b,lb,ub,x0)同上，并设置初值x0。 x = quadprog(H,f,A,b,lb,ub,x0,options)根据options参数指定的优化参数进行最小化。

[x,fval] = quadprog(...)返回解x处的目标函数值fval = 0.5*x'*H*x + f'*x。 [x,fval,exitflag] = quadprog(...)返回exitflag参数，描述计算的退出条件。 [x,fval,exitflag,output] = quadprog(...)返回包含优化信息的结构输出output。

[x,fval,exitflag,output,lambda] = quadprog(...)返回解x处包含拉格朗日乘子的lambda参数。 变量：

各变量的意义同前。 注意： 1．

一般地，如果问题不是严格凸性的，用quadprog函数得到的可能是局部最优解。 2．

如果用Aeq和Beq明确地指定等式约束，而不是用lb和ub指定，则可以得到更好的数值解。 3．

若x的组分没有上限或下限，则quadprog函数希望将对应的组分设置为Inf（对于上限）或-Inf（对于下限），而不是强制性地给予上限一个很大的数或给予下限一个很小的负数。 4．

对于大型优化问题，若没有提供初值x0，或x0不是严格可行，则quadprog函数会选择一个新的初始可行点。 5．

若为等式约束，且quadprog函数发现负曲度(negative curvature)，则优化过程终止，exitflag的值等于-1。 算法：

大型优化算法

当优化问题只有上界和下界，而没有线性不等式或等式约束，则缺省算法为大型算法。或者，如果优化问题中只有线性等式，而没有上界和下界或线性不等式时，缺省算法也是大型算法。

本法是基于内部映射牛顿法(interior-reflective Newton method)的子空间置信域法(subspace trust-region)。该法的具体算法请参见文献[2]。该法的每一次迭代都与用PCG法求解大型线性系统得到的近似解有关。 中型优化算法

quadprog函数使用活动集法，它也是一种投影法，首先通过求解线性规划问题

来获得初始可行解。 诊断： 1.

大型优化问题

大型优化问题不允许约束上限和下限相等，如若lb(2)==ub(2)，则给出以下出错信息：

Equal upper and lower bounds not permitted in this large-scale method.Use equality constraints and the medium-scale method instead.

若优化模型中只有等式约束，仍然可以使用大型算法；如果模型中既有等式约束又有边界约束，则必须使用中型方法。 2．中型优化问题

当解不可行时，quadprog函数给出以下警告：

Warning: The constraints are overly stringent;there is no feasible solution.

这里，quadprog函数生成使约束矛盾最坏程度最小的结果。

当等式约束不连续时，给出下面的警告信息：

Warning: The equality constraints are overly stringent;there is no feasible solution.

当Hessian矩阵为负半定时，生成无边界解，给出下面的警告信息：

Warning: The solution is unbounded and at infinity;the constraints are not restrictive enough.

这里，quadprog函数返回满足约束条件的x值。 局限性： 1．

此时，显示水平只能选择'off'和'final'，迭代参数'iter'不可用。 2．

当问题不定或负定时，常常无解（此时exitflag参数给出一个负值，表示优化过程不收敛）。若正定解存在，则quadprog函数可能只给出局部极小值，因为问题可能时非凸的。 3．

对于大型问题，不能依靠线性等式，因为Aeq必须是行满秩的，即Aeq的行数必须不多于列数。若不满足要求，必须调用中型算法进行计算。

9.2.4.3 应用实例 [例一]

求解下面的最优化问题：

目标函数 约束条件

首先，目标函数可以写成下面的矩阵形式： 其中

输入下列系数矩阵： H = [1 -1; -1 2] f = [-2; -6]

A = [1 1; -1 2; 2 1] b = [2; 2; 3] lb = zeros(2,1)

然后调用二次规划函数quadratic：

[x,fval,exitflag,output,lambda] = quadprog(H,f,A,b,[],[],lb) 得问题的解 x =

0.6667

1.3333 fval =

-8.2222 exitflag = 1

output =

iterations: 3

algorithm: 'medium-scale: active-set'

firstorderopt: []

cgiterations: [] lambda.ineqlin ans =

3.1111

0.4444

0

lambda.lower ans = 0 0

磁盘中本问题的M文件为opt24_1.m。

9.2.5 有约束最小化 9.2.5.1 基本数学原理 在有约束最优化问题中，通常要将该问题转换为更简单的子问题，这些子问题可以求解并作为迭代过程的基础。早期的方法通常是通过构造惩罚函数等来将有约束的最优化问题转换为无约束最优化问题进行求解。现在，这些方法已经被更有效的基于K-T（Kuhn-Tucker）方程解的方法所取代。K-T方程是有约束最优化问题求解的必要条件。假设有所谓的Convex规划问题，f(x)和Gi(x)，i=1,2,…,m为Convex函数，则K-T方程对于求得全局极小点是必要的，也是充分的。 对于规划问题

其中，x是设计参数向量，（x∈Rn），f(x)为目标函数，返回标量值，向量函数G(x)返回等式约束和不等式约束在x处的值。 它的K-T方程可表达为：

（*）

i=1,?,m

i=me+1,?,m

其中第一行描述了目标函数和约束条件在解处梯度的取消。由于梯度取消，需要用拉格朗日乘子λi(i=1,2,…,m)来平衡目标函数与约束梯度间大小的差异。

K-T方程的解形成了许多非线性规划算法的基础，这些算法直接计算拉格朗日乘子，通过用拟牛顿法更新过程，给K-T方程积累二阶信息，可以保证有约束拟牛顿法的超线性收敛。这些方法称为序列二次规划法（SQP），因为在每一次主要的迭代中都求解一次二次规划问题。

对于给定的规划问题，序列二次规划（SQP）的主要的思路是形成基于拉格朗日函数二次近似的二次规划子问题，即

这里，通过假设约束条件为不等式约束来使（*）式得到了简化，通过非线性有约束问题线性化来获得二次规划子问题。

二次规划子问题可表达为

i=1,…,me

i=me+1,…,m

该子问题可以用任意一种二次规划算法求解，求得的解可以用来形成新的迭代公式xk+1=xk+αkdk。

用SQP法求解非线性有约束问题时的迭代次数常比用解无约束问题时的少，因为在搜索区域内，SQP方法可以获得最佳的搜索方向和步长信息。 给Rosenbrock函数添加非线性不等式约束g(x)

经过96次迭代得到问题的解：x=[0.9072,0.8288]，初值为x=[-1.9,2]，无约束问题则需要140次迭代。图9-3是搜索路径图。

图9-3 SQP法的搜索路径

MATLAB中SQP法的实现分三步，即

●

拉格朗日函数Hessian矩阵的更新；

●

二次规划问题求解；

●

一维搜索和目标函数的计算 (一)、Hessian矩阵的更新

在每一次主要迭代过程中，都用BFGS法计算拉格朗日函数的Hessian矩阵的拟牛顿近似矩阵。更新公式为： 其中：

上式中，λi(i=1,2,…,m)为拉格朗日乘子的估计。 （二）、二次规划求解

SQP法的每一次主要迭代过程中都要求一次二次规划问题，形式如下：

i=1,…,me

i=me+1,…,m

求解过程分两步，第一步涉及可行点（若存在）的计算，第二步为可行点至解的迭代序列。

在第一步中，需要有可行点作为初值，若当前点不可行，则通过求解下列线性规划问题可以得到一个可行点：

i=1,…,me

i=me+1,…,m

其中，Ai为矩阵A的第i行。

（三）、一维搜索和目标函数的计算

二次规划子问题的解生成一个向量dk，它形成一个新的迭代公式： αk为步长参数。

目标函数的形式如下： 其中，

i=1,…,m

9.2.5.2 相关函数介绍 fmincon函数

功能：求多变量有约束非线性函数的最小值。

数学模型：

其中，x, b, beq, lb,和ub为向量， A 和 Aeq 为矩阵， c(x) 和 ceq(x)为函数，返回标量。f(x), c(x), 和 ceq(x)可以是非线性函数。 语法：

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 = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options,P1,P2, ...) [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(...) 描述：

fmincon 求多变量有约束非线性函数的最小值。该函数常用于有约束非线性优化问题。

x = fmincon(fun,x0,A,b) 给定初值x0，求解fun函数的最小值x。fun函数的约束条件为A*x <= b，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的下界lb和上界ub，使得总是有lb <= x <= ub。若无等式存在，则令Aeq=[]、beq=[]。 x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon) 在上面的基础上，在

nonlcon参数中提供非线性不等式c(x)或等式ceq(x)。 fmincon函数要求c(x) <= 0且ceq(x) = 0。当无边界存在时，令lb=[]和（或）ub=[]。

x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options) 用optiions参数指定的参数进行最小化。

x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options,P1,P2,...) 将问题参数P1, P2等直接传递给函数fun和nonlin。若不需要这些变量，则传递空矩阵到A, b, Aeq, beq, lb, ub, nonlcon和 options。 [x,fval] = fmincon(...) 返回解x处的目标函数值。 [x,fval,exitflag] = fmincon(...) 返回exitflag参数，描述函数计算的退出条件。

[x,fval,exitflag,output] = fmincon(...) 返回包含优化信息的输出参数output。

[x,fval,exitflag,output,lambda] = fmincon(...) 返回解x处包含拉格朗日乘子的lambda参数。

[x,fval,exitflag,output,lambda,grad] = fmincon(...) 返回解x处fun函数的梯度。

[x,fval,exitflag,output,lambda,grad,hessian] = fmincon(...) 返回解x处fun函数的Hessian矩阵。 变量：

nonlcon参数

该参数计算非线性不等式约束c(x)<=0 和非线性等式约束ceq(x)=0。 nonlcon 参数是一个包含函数名的字符串。该函数可以是M文件、内部文件或MEX文件。它要求输入一个向量x，返回两个变量—解x处的非线性不等式向量c和非线性等式向量ceq。例如，若nonlcon='mycon'，则M文件mycon.m具有下面的形式： function [c,ceq] = mycon(x) c = ...

% 计算x处的非线性不等式。 ceq = ...

% 计算x处的非线性等式。 若还计算了约束的梯度，即

options = optimset('GradConstr','on')

则nonlcon函数必须在第三个和第四个输出变量中返回c(x)的梯度GC和ceq(x)的梯度Gceq。当被调用的nonlcon函数只需要两个输出变量（此时优化算法只需要c和ceq的值，而不需要GC和GCeq）时，可以通过查看nargout的值来避免计算GC和Gceq的值。

function [c,ceq,GC,GCeq] = mycon(x)

c = ...

% 解x处的非线性不等式。

ceq = ...

% 解x处的非线性等式。

if nargout > 2

% 被调用的nonlcon函数，要求有4个输出变量。

GC = ...

% 不等式的梯度。

GCeq = ...

% 等式的梯度。 end

若nonlcon函数返回m元素的向量c和长度为n的x，则c(x)的梯度GC是一个n*m的矩阵，其中GC(i,j)是c(j)对x(i)的偏导数。同样，若ceq是一个p元素的向量，则ceq(x)的梯度Gceq是一个n*p的矩阵，其中Gceq(i,j)是ceq(j)对x(i)的偏导数。

其它参数意义同前。 注意：

大型优化问题：

1． 使用大型算法，必须在fun函数中提供梯度信息(options.GradObj设置为'on')。如果没有梯度信息，则给出警告信息。

fmincon函数允许g(x)为一近似梯度，但使用真正的梯度将使优化过程更具稳健性。

2．当对矩阵的二阶导数（即Hessian矩阵）进行计算以后，用该函数求解大型问题将更有效。但不需要求得真正的Hessian矩阵，如果能提供Hessian矩阵的稀疏结构的信息（用options参数的HessPattern属性），则fmincon函数可以算得Hessian矩阵的稀疏有限差分近似。

3．若x0不是严格可行的，则fmincon函数选择一个新的严格可行初始点。 4．若x的某些元素没有上界或下界，则fmincon函数更希望对应的元素设置为Inf（对于上界）或-Inf（对于下界），而不希望强制性地给上界赋一个很大的值，给下界赋一个很小的负值。

5．线性约束最小化课题中也有几个问题需要注意： ??

Aeq矩阵中若存在密集列或近密集列（A dense (或fairly dense)column），会导致满秩并使计算费时。 ??

fmincon函数剔除Aeq中线性相关的行。此过程需要进行反复的因子分解，因此，如果相关行很多的话，计算将是一件很费时的事情。 ??

每一次迭代都要用下式进行稀疏最小二乘求解

其中RT为前提条件的乔累斯基因子。 中型优化问题 1．如果用Aeq和beq清楚地提供等式约束，将比用lb和ub获得更好的数值解。 2．在二次子问题中，若有等式约束并且因等式（dependent equalities）被发现和剔除的话，将在过程标题中显示'dependent'(当output参数要求使用options.Display = 'iter')。只有在等式连续的情况下，因等式才会被剔除。若等式系统不连续，则子问题将不可行并在过程标题中打印'infeasible'信息。 算法：

大型优化算法 缺省时，若提供了函数的梯度信息，并且只有上下界存在或只有线性等式约束存在，则fmincon函数将选择大型算法。本法是基于内部映射牛顿法(interior-reflective Newton method)的子空间置信域法(subspace

trust-region)。该法的具体算法请参见文献[2]。该法的每一次迭代都与用PCG法求解大型线性系统得到的近似解有关。 中型优化算法

fmincon函数使用序列二次规划法(SQP)。本法中，在每一步迭代中求解二次规划子问题，并用BFGS法更新拉格朗日Hessian矩阵。参见文献[3]、[6]。 该算法使用与文献[1]、[2]和[3]中提供的相似的目标函数进行一维搜索。二次子问题使用文献[4]描述的活动集方案进行求解。 诊断：

大型优化问题

求大型优化问题的代码中不允许上限和下限相等，即不能有lb(2)==ub(2),否则给出下面的出错信息：

Equal upper and lower bounds not permitted in this large-scale

method.

Use equality constraints and the medium-scale method instead.

若只有等式约束，可以仍然使用大型算法。当既有等式约束又有边界约束时，使用中型算法。 局限： 1．

目标函数和约束函数都必须是连续的，否则可能会给出局部最优解。 2．

当问题不可行时，fmincon函数将试图使最大约束值最小化。 3．

目标函数和约束函数都必须是实数。 4．

对于大型优化问题，使用大型优化算法时，用户必须在fun函数中提供梯度(options参数的GradObj属性必须设置为'on') , 并且只可以指定上界和下界约束，或者只有线性约束必须存在，Aeq的行数不能多于列数。 5． 现在，如果在fun函数中提供了解析梯度，选项参数DerivativeCheck不能与大

解，此时必须引进非劣解的概念（非劣解又称为有效解或帕累托解）。 定义：若x*(x*∈Ω)的邻域内不存在Δx，使得（x*+Δx）∈Ω，且 i=1,…,m

对于某些j

则称为x*非劣解。

多目标规划有许多解法，下面列出常用的几种： 1. 权和法

该法将多目标向量问题转化为所有目标的加权求和的标量问题，即

加权因子的选取方法很多，有专家打分法、α方法、容限法和加权因子分解法等。

该问题可以用标准的无约束最优化算法进行求解。 2.

ε约束法

ε约束法克服了权和法的某些凸性问题。它对目标函数向量中的主要目标Fp进行最小化，将其它目标用不等式约束的形式写出： sub.

i=1,…,m i≠p 3.

目标达到法

目标函数系列为F(x)={F1(x), F2(x),…, Fm(x)}，对应地有其目标值系列 。允许目标函数有正负偏差，偏差的大小由加权系数向量W={W1,W2,…,Wm}控制，于是目标达到问题可以表达为标准的最优化问题: sub.

i=1,…,m

指定目标{ }，定义目标点P。权重向量定义从P到可行域空间Λ（γ）的搜索方向，在优化过程中，γ的变化改变可行域的大小，约束边界变为唯一解点F1s、F2s。

4.目标达到法的改进

目标达到法的一个好处是可以将多目标最优化问题转化为非线性规划问题，但是，在序列二次规划（SQP）过程中，一维搜索的目标函数选择不是一件容易的事情，因为在很多情况下，很难决定是使目标函数变大好还是使它变小好。这导致许多目标函数创建过程的提出。可以通过将目标达到问题变为最大最小化问题来获得更合适的目标函数。 其中

i=1,…,m

9.2.6.2 相关函数介绍 fgoalattain函数

功能：求解多目标达到问题

数学模型： 其中x, weight, goal, b, beq, lb和ub 为向量, A和Aeq为矩阵, c(x), ceq(x)和F(x)为函数，返回向量。F(x), c(x)和ceq(x)可以是非线性函数。 语法：

x = fgoalattain(fun,x0,goal,weight) x = fgoalattain(fun,x0,goal,weight,A,b)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon) x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,...

lb,ub,nonlcon,options)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,...

lb,ub,nonlcon,options,P1,P2,...) [x,fval] = fgoalattain(...)

[x,fval,attainfactor] = fgoalattain(...)

[x,fval,attainfactor,exitflag] = fgoalattain(...)

[x,fval,attainfactor,exitflag,output] = fgoalattain(...)

[x,fval,attainfactor,exitflag,output,lambda] = fgoalattain(...) 描述：

fgoalattain求解多目标达到问题。

x = fgoalattain(fun,x0,goal,weight)试图通过变化x来使目标函数fun达到goal指定的目标。初值为x0,weight参数指定权重。

x = fgoalattain(fun,x0,goal,weight,A,b)求解目标达到问题，约束条件为线性不等式A*x <= b。

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq)求解目标达到问题，除提供上面的线性不等式以外，还提供线性等式Aeq*x = beq 。当没有不等式存在时，设置A=[]、b=[]。

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub)为设计变量x定义下界lb和上界ub集合，这样始终有lb <= x <= ub。 x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon) 将目标达到问题归结为nonlcon参数定义的非线性不等式c(x)或非线性等式ceq(x)。fgoalattain函数优化的约束条件为(x) <= 0和ceq(x) = 0。若不存在边界，设置lb=[]和（或）ub=[]。

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon,... options)用options中设置的优化参数进行最小化。

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon,...

options,P1,P2,...) 将问题参数P1, P2等直接传递给函数fun和nonlcon。如果不需要参数A, b, Aeq, beq, lb, ub, nonlcon和options，将它们设置为空矩阵。

[x,fval] = fgoalattain(...) 返回解x处的目标函数值。

[x,fval,attainfactor] = fgoalattain(...) 返回解x处的目标达到因子。 [x,fval,attainfactor,exitflag] = fgoalattain(...)返回exitflag参数，描

述计算的退出条件。

[x,fval,attainfactor,exitflag,output] = fgoalattain(...)返回包含优化信息的输出参数output。

[x,fval,attainfactor,exitflag,output,lambda] = fgoalattain(...) 返回包含拉格朗日乘子的lambda参数。 变量： goal变量

目标希望达到的向量值。向量的长度与fun函数返回的目标数F相等。

fgoalattain函数试图通过最小化向量F中的值来达到goal参数给定的目标。 nonlcon参数：

该函数计算非线性不等式约束c(x) <=0和非线性等式约束ceq(x)=0。nonlcon函数是一个包含函数名的字符串，该函数可以是M文件、内部函数或MEX文件。nonlcon函数需要输入向量x，返回两个变量—x处的非线性不等式向量c和x处的非线性等式向量ceq。例如，若nonlcon='mycon'，则M文件的形式如下： function [c,ceq] = mycon(x) c = ...

% 计算x处的非线性不等式。 ceq = ...

% 计算x处的非线性等式。

若约束函数的梯度可以计算，且options.GradConstr设为'on',即

options = optimset('GradConstr','on')

则函数nonlcon也必须在第三个和第四个输出变量中输出c(x)的梯度GC和

ceq(x)的梯度GCeq。注意，可以通过核对nargout参数来避免计算GC和GCeq。 function [c,ceq,GC,GCeq] = mycon(x) c = ...

% x处的非线性不等式。 ceq = ...

% x处的非线性等式。 if nargout > 2

% 被调用的nonlcon函数，有4个输出。

GC = ...

% 不等式的梯度。

GCeq = ...

% 等式的梯度。 end

若nonlcon函数返回m元素的向量c和长度为n的x，则c(x)的梯度GC是一个n*m的矩阵，其中GC(i,j)是c(j)对x(i)的偏导数。同样，若ceq是一个p元素的向量，则ceq(x)的梯度Gceq是一个n*p的矩阵，其中Gceq(i,j)是ceq(j)对x(i)的偏导数。

options变量

优化参数选项。你可以用optimset函数设置或改变这些参数的值。 DerivativeCheck – 比较用户提供的导数（目标函数或约束函数的梯度）和有限差分导数。

Diagnostics – 打印将要最小化或求解的函数的诊断信息。 DiffMaxChange – 变量中有限差分梯度的最大变化。 DiffMinChange - 变量中有限差分梯度的最小变化。

Display – 显示水平。设置为'off'时不显示输出；设置为'iter'时显示每一次迭代的输出；设置为'final'时只显示最终结果。

GoalExactAchieve – 使得目标个数刚好达到，不多也不少。 GradConstr – 用户定义的约束函数的梯度。 GradObj – 用户定义的目标函数的梯度。使用大型方法时必须使用梯度，对于中型方法则是可选项。

MaxFunEvals – 函数评价的允许最大次数。 MaxIter – 函数迭代的允许最大次数。 MeritFunction – 如果设为'multiobj'，则使用目标达到或最大最小化目标函数的方法。若设置为'singleobj'，则使用fmincon 函数计算目标函数。 TolCon – 约束矛盾的终止容限。 TolFun – 函数值处的终止容限。 TolX – x处的终止容限。

weight变量

为权重向量，可以控制低于或超过fgoalattain函数指定目标的相对程度。当goal的值都是非零值时，为了保证活动对象超过或低于的比例相当，将权重函数设置为abs(goal) (活动对象为阻止解处目标改善的对象集合)。

注意: 1.

当目标值中的任意一个为零时，设置weight=abs(goal)将导致目标约束看起来更象硬约束，而不象目标约束。 2.

当加权函数weight为正时，fgoalattain函数试图使对象小于目标值。为了使目标函数大于目标值，将权重weight设置为负。为了使目标函数尽可能地接近目标值，使用GoalsExactAchieve参数，将fun函数返回的第一个元素作为目标。

attainfactor变量

attainfactor变量是超过或低于目标的个数。若attainfactor为负，则目标已经溢出；若attainfactor为正，则目标个数还未达到。

其它参数意义同前。 注意：

当特征值为复数时，本问题不连续，这也说明了为什么收敛速度很慢。尽管原始方法假设函数是连续的，该法仍然可以向解的方向前进，因为在解的位置上，没有发生不连续的现象。当对象和目标为复数时，fgoalattain函数将试图得到最小二乘意义上的目标。

算法：

多目标优化同时涉及到一系列对象。fgoalattain函数求解该问题的基本算法是目标达到法。该法为目标函数建立起目标值。多目标优化的具体算法在前面进行了更详细的介绍。

本次实现过程中，使用了松弛变量γ作为模糊变量同时最小化目标向量F(x)；goal参数是一系列目标达到值。在进行优化之前，通常不知道对象是否会达到目标。使用权向量weight可以控制是没有达到还是溢出。

fgoalattain函数使用序列二次规划法(SQP)，前面已经进行了比较多的介绍。算法中对于一维搜索和Hessian矩阵进行了修改。在一维搜索中，将精确目标函数和文献[2]、[3]中的目标函数一起使用。当有一个目标函数不再发生改善时，一维搜索终止。修改的Hessian矩阵借助于本问题的结构，也被采用。详细内容可参见文献[5]、[6]。

attainfactor参数包含解处的γ值。γ取负值时表示目标溢出。 局限性：

目标函数必须是连续的。fgoalattain函数将只给出局部最优解。 参见：

fmincon, fminimax, optimset 9.2.6.3 应用实例

[例二] 某化工厂拟生产两种新产品A和B，其生产设备费用分别为：A，2万元/吨；B，5万元/吨。这两种产品均将造成环境污染，设由公害所造成的损失可折算为：A，4万元/吨；B，1万元/吨。由于条件限制，工厂生产产品A和B的最大生产能力各为每月5吨和6吨，而市场需要这两种产品的总量每月不少于7吨。试问工厂如何安排生产计划，在满足市场需要的前提下，使设备投资和公害损失均达最小。该工厂决策认为，这两个目标中环境污染应优先考虑，设备投资的目标值为20万元，公害损失的目标为12万元。

设工厂每月生产产品A为x1吨，B为x2吨，设备投资费为f1(x)，公害损失费为f2(x)，则这个问题可表达为多目标优化问题：

首先需要编写目标函数的M文件opt26_2o.m，返回目标计算值

function f=myfun(x)

f(1)=2*x(1)+5*x(2);

f(2)=4*x(1)+x(2);

给定目标，权重按目标比例确定，给出初值

goal=[20 12];

weight=[20 12];

x0=[2 5];

给出约束条件的系数

本文来源：https://www.bwwdw.com/article/i2zt.html