Main Content

编写标量目标函数

函数文件

标量目标函数文件接受一个输入,例如 x,并返回一个实数标量输出,例如 f。输入 x 可以是标量、向量或矩阵。一个函数文件可以返回多个输出(请参阅包括梯度和黑塞函数)。

例如,假设您的目标是包含三个变量(x、y 和 z)的函数:

f(x) = 3*(x – y)4 + 4*(x + z)2 / (1 + x2 + y2 + z2) + cosh(x – 1) + tanh(y + z)。

  1. 将此函数编写为接受向量 xin = [x;y;z] 并返回 f 的文件:

    function f = myObjective(xin)
    f = 3*(xin(1)-xin(2))^4 + 4*(xin(1)+xin(3))^2/(1+norm(xin)^2) ...
        + cosh(xin(1)-1) + tanh(xin(2)+xin(3));
  2. 将其以名为 myObjective.m 的文件保存到 MATLAB® 路径上的文件夹中。

  3. 检查函数是否能够正确进行计算:

    myObjective([1;2;3])
    
    ans =
        9.2666

有关如何包含额外参数的信息,请参阅传递额外参数。有关函数文件的更复杂示例,请参阅Minimization with Gradient and Hessian Sparsity PatternMinimization with Bound Constraints and Banded Preconditioner

局部函数和嵌套函数

函数可以作为局部函数嵌套函数存在于其他文件中。使用局部函数或嵌套函数可以减少保存的不同文件的数量。使用嵌套函数还允许您访问额外的参数,如嵌套函数中所示。

例如,假设您要在非线性约束中所述的 ellipseparabola.m 约束条件下对函数文件中所述的 myObjective.m 目标函数进行最小化。您不必编写两个文件 myObjective.mellipseparabola.m,而只需编写一个将这两个函数作为局部函数包含在其中的文件:

function [x fval] = callObjConstr(x0,options)
% Using a local function for just one file

if nargin < 2
    options = optimoptions('fmincon','Algorithm','interior-point');
end

[x fval] = fmincon(@myObjective,x0,[],[],[],[],[],[], ...
    @ellipseparabola,options);

function f = myObjective(xin)
f = 3*(xin(1)-xin(2))^4 + 4*(xin(1)+xin(3))^2/(1+sum(xin.^2)) ...
    + cosh(xin(1)-1) + tanh(xin(2)+xin(3));

function [c,ceq] = ellipseparabola(x)
c(1) = (x(1)^2)/9 + (x(2)^2)/4 - 1;
c(2) = x(1)^2 - x(2) - 1;
ceq = [];

从点 [1;1;1] 开始求解有约束最小化问题:

[x fval] = callObjConstr(ones(3,1))

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is 
non-decreasing in feasible directions, to within the default 
value of the function tolerance, and constraints are satisfied 
to within the default value of the constraint tolerance.

x =
    1.1835
    0.8345
   -1.6439

fval =
    0.5383

匿名函数目标

使用匿名函数编写简单的目标函数。有关匿名函数的详细信息,请参阅什么是匿名函数?。罗森布罗克函数非常简单,可以写成匿名函数:

anonrosen = @(x)(100*(x(2) - x(1)^2)^2 + (1-x(1))^2);
检查 anonrosen[-1 2] 处的计算是否正确:
anonrosen([-1 2])

ans =
   104
使用 fminunc 最小化 anonrosen 会生成以下结果:
options = optimoptions(@fminunc,'Algorithm','quasi-newton');
[x fval] = fminunc(anonrosen,[-1;2],options)

Local minimum found.

Optimization completed because the size of the gradient
is less than the default value of the function tolerance.

x =
    1.0000
    1.0000

fval =
  1.2266e-10

包括梯度和黑塞函数

为求解器提供导数

对于 fminconfminunc,您可以在目标函数中包含梯度。一般情况下,当包含梯度时,求解器更稳健,速度也会有所提升。请参阅包含导数的好处。要同时包括二阶导数(黑塞矩阵),请参阅包含黑塞函数

下表显示哪些算法可以使用梯度和黑塞函数。

求解器算法梯度黑塞矩阵
fminconactive-set可选
interior-point可选可选(请参阅适用于 fmincon 内点算法的黑塞函数
sqp可选
trust-region-reflective必需可选(请参阅适用于 fminunc 信赖域或 fmincon 信赖域反射算法的黑塞函数
fminuncquasi-newton可选
trust-region必需可选(请参阅适用于 fminunc 信赖域或 fmincon 信赖域反射算法的黑塞函数

如何包含梯度

  1. 编写返回以下内容的代码:

    • 目标函数(标量)作为第一个输出

    • 梯度(向量)作为第二个输出

  2. 使用 optimoptionsSpecifyObjectiveGradient 选项设置为 true。适用时,也将 SpecifyConstraintGradient 选项设置为 true

  3. (可选)检查梯度函数是否与有限差分逼近匹配。请参阅Checking Validity of Gradients or Jacobians

提示

为了获得最大的灵活性,请编写条件化代码。条件化意味着函数输出的数目可以变化,如以下示例中所示。条件化代码不会因 SpecifyObjectiveGradient 选项的值而出错。非条件化代码要求您适当地设置选项。

例如,假设有以下罗森布罗克函数

f(x)=100(x2x12)2+(1x1)2,

使用优化实时编辑器任务或求解器的有约束非线性问题中对该函数进行了说明和绘图。f(x) 的梯度为

f(x)=[400(x2x12)x12(1x1)200(x2x12)],

rosentwo 是条件化函数,可返回求解器所需的任何内容:

function [f,g] = rosentwo(x)
% Calculate objective f
f = 100*(x(2) - x(1)^2)^2 + (1-x(1))^2;

if nargout > 1 % gradient required
    g = [-400*(x(2)-x(1)^2)*x(1)-2*(1-x(1));
        200*(x(2)-x(1)^2)];
    
end

nargout 检查调用函数指定的参量个数。请参阅查找函数参量的数量

为无约束优化设计的 fminunc 求解器允许您最小化罗森布罗克函数。通过设置 optionsfminunc 使用梯度和黑塞函数:

options = optimoptions(@fminunc,'Algorithm','trust-region',...
    'SpecifyObjectiveGradient',true);

[-1;2] 开始运行 fminunc

[x fval] = fminunc(@rosentwo,[-1;2],options)
Local minimum found.

Optimization completed because the size of the gradient
is less than the default value of the function tolerance.

x =
    1.0000
    1.0000

fval =
  1.9886e-17

如果您有 Symbolic Math Toolbox™ 许可证,您可以自动计算梯度和黑塞矩阵,如Calculate Gradients and Hessians Using Symbolic Math Toolbox中所述。

包含黑塞函数

您可以使用 fmincon 'trust-region-reflective''interior-point' 算法以及 fminunc 'trust-region' 算法来包含二阶导数。根据信息的类型和算法,可通过几种方法来包括黑塞函数信息。

您还必须包含梯度(将 SpecifyObjectiveGradient 设置为 true,如果适用,还必须将 SpecifyConstraintGradient 设置为 true)以便包含黑塞函数。

适用于 fminunc 信赖域或 fmincon 信赖域反射算法的黑塞函数.  这些算法要么没有约束,要么只有边界或线性等式约束。因此,黑塞矩阵是由目标函数的二阶导数组成的矩阵。

包含黑塞函数作为目标函数的第三个输出。例如,罗森布罗克函数的黑塞函数 H(x) 如下(请参阅如何包含梯度

H(x)=[1200x12400x2+2400x1400x1200].

在目标中包含此黑塞函数:

function [f, g, H] = rosenboth(x)
% Calculate objective f
f = 100*(x(2) - x(1)^2)^2 + (1-x(1))^2;

if nargout > 1 % gradient required
    g = [-400*(x(2)-x(1)^2)*x(1)-2*(1-x(1));
        200*(x(2)-x(1)^2)];
    
    if nargout > 2 % Hessian required
        H = [1200*x(1)^2-400*x(2)+2, -400*x(1);
            -400*x(1), 200];  
    end

end

HessianFcn 设置为 'objective'。例如,

options = optimoptions('fminunc','Algorithm','trust-region',...
    'SpecifyObjectiveGradient',true,'HessianFcn','objective');

适用于 fmincon 内点算法的黑塞函数.  该黑塞函数是拉格朗日函数的黑塞函数,其中 L(x,λ) 是

L(x,λ)=f(x)+λg,igi(x)+λh,ihi(x).

g 和 h 是向量函数,分别表示所有不等式和等式约束(指有界、线性和非线性约束),因此最小化问题的公式为

minxf(x) subject to g(x)0, h(x)=0.

有关详细信息,请参阅受约束的最优性理论。拉格朗日方程的黑塞函数公式为

xx2L(x,λ)=2f(x)+λg,i2gi(x)+λh,i2hi(x).(1)

要包含黑塞函数,请用以下语法编写函数

hessian = hessianfcn(x,lambda)

hessian 是 n×n 矩阵(稀疏或稠密),其中 n 是变量的数目。如果 hessian 很大并且非零项相对较少,请将 hessian 表示为稀疏矩阵,以节省运行时间和内存。lambda 是具有与非线性约束相关联的拉格朗日乘数向量的结构体:

lambda.ineqnonlin
lambda.eqnonlin

fmincon 计算结构体 lambda,并将其传递给您的黑塞函数。hessianfcn 必须计算公式 1 中的总和。通过设置以下选项表明您提供了黑塞函数:

options = optimoptions('fmincon','Algorithm','interior-point',...
    'SpecifyObjectiveGradient',true,'SpecifyConstraintGradient',true,...
    'HessianFcn',@hessianfcn);

例如,要针对罗森布罗克函数包含黑塞函数(约束为单位圆盘 x12+x221),请注意约束函数 g(x)=x12+x2210 具有梯度和二阶导数矩阵

g(x)=[2x12x2]Hg(x)=[2002].

将黑塞函数写成

function Hout = hessianfcn(x,lambda)
% Hessian of objective
H = [1200*x(1)^2-400*x(2)+2, -400*x(1);
            -400*x(1), 200];
% Hessian of nonlinear inequality constraint
Hg = 2*eye(2);
Hout = H + lambda.ineqnonlin*Hg;

hessianfcn 保存到 MATLAB 路径。为了完成此示例,包含梯度的约束函数为

function [c,ceq,gc,gceq] = unitdisk2(x)
c = x(1)^2 + x(2)^2 - 1;
ceq = [ ];

if nargout > 2
    gc = [2*x(1);2*x(2)];
    gceq = [];
end

求解包括梯度和黑塞函数的问题。

fun = @rosenboth;
nonlcon = @unitdisk2;
x0 = [-1;2];
options = optimoptions('fmincon','Algorithm','interior-point',...
    'SpecifyObjectiveGradient',true,'SpecifyConstraintGradient',true,...
    'HessianFcn',@hessianfcn);
[x,fval,exitflag,output] = fmincon(fun,x0,[],[],[],[],[],[],@unitdisk2,options);

有关使用内点黑塞函数的其他示例,请参阅使用解析黑塞函数的 fmincon 内点算法Calculate Gradients and Hessians Using Symbolic Math Toolbox

黑塞矩阵乘法函数.  fmincon interior-pointtrust-region-reflective 算法都允许您提供黑塞函数乘法函数,而无需提供完整的黑塞函数。此函数给出黑塞函数乘以向量的乘积结果,而不直接计算黑塞函数。这可以节省内存。SubproblemAlgorithm 选项必须为 'cg',黑塞矩阵乘法函数才能发挥作用;这是 trust-region-reflective 的默认值。

这两种算法的语法不同。

  • 对于 interior-point 算法,语法如下

    W = HessMultFcn(x,lambda,v);

    结果 W 应为 H*v 的积,其中 H 是在 x 处拉格朗日方程的黑塞函数(请参阅公式 1),lambda 是拉格朗日乘数(由 fmincon 计算),v 是大小为 n×1 的向量。按如下所示设置选项:

    options = optimoptions('fmincon','Algorithm','interior-point','SpecifyObjectiveGradient',true,... 
        'SpecifyConstraintGradient',true,'SubproblemAlgorithm','cg','HessianMultiplyFcn',@HessMultFcn);

    提供函数 HessMultFcn,该函数返回 n×1 向量,其中 n 是 x 的维数。HessianMultiplyFcn 选项使您能够传递黑塞函数乘以向量的结果,而无需计算黑塞函数。

  • trust-region-reflective 算法不涉及 lambda

    W = HessMultFcn(H,v);

    结果 W = H*vfmincon 传递 H 作为在目标函数的第三个输出中返回的值(请参阅适用于 fminunc 信赖域或 fmincon 信赖域反射算法的黑塞函数)。fmincon 还传递 v,即包含 n 个行的向量或矩阵。v 中的列数可能会有所不同,因此请写入 HessMultFcn 以接受任意数量的列。H 不一定是黑塞函数;它可以是任何让您能够计算 W = H*v 的项。

    按如下所示设置选项:

    options = optimoptions('fmincon','Algorithm','trust-region-reflective',... 
        'SpecifyObjectiveGradient',true,'HessianMultiplyFcn',@HessMultFcn);

    有关将黑塞函数乘法函数和 trust-region-reflective 算法结合使用的示例,请参阅Minimization with Dense Structured Hessian, Linear Equalities

包含导数的好处

如果您不提供梯度,求解器会通过有限差分来估计梯度。如果提供梯度,则求解器不需要执行这种有限差分估计,从而可以节省时间且结果更准确,但对于复杂导数来说,有限差分估计可能更快。此外,求解器使用逼近的黑塞函数,它可能与真实的黑塞函数相差甚远。提供黑塞函数可以用更少的迭代次数来生成解。有关示例,请参阅Calculate Gradients and Hessians Using Symbolic Math Toolbox 的末尾部分。

对于约束问题,提供梯度还有另一个优势。求解器可以到达一个点 x,这说明 x 是可行点,但对于此处的 xx 附近的有限差分总是导致不可行点。进一步假设在不可行点处的目标函数返回复数输出、InfNaN 或错误。在这种情况下,求解器可能会失败或过早停止。提供梯度可允许求解器继续。要获得此优势,您可能还需要包括非线性约束函数的梯度,并将 SpecifyConstraintGradient 选项设置为 true。请参阅非线性约束

为内点 fmincon 选择输入黑塞函数逼近方法

fmincon interior-point 算法提供了许多选项用来选择输入黑塞函数逼近方法。有关语法详细信息,请参阅作为输入的黑塞矩阵。以下是选项以及对其相对特征的估计。

黑塞矩阵相对内存使用量相对效率
'bfgs'(默认值)高(适用于大型问题)
'lbfgs'低到中等中等
'fin-diff-grads'中等
'HessianMultiplyFcn'低(取决于您的代码)中等
'HessianFcn'不确定(取决于您的代码)高(取决于您的代码)

请使用默认的 'bfgs'黑塞函数,除非您

'lbfgs' 只有中等效率的原因是需要进行双倍计算。它采用计算成本相对高昂的谢尔曼-莫里森更新。由于 'lbfgs' 有内存使用限制,得到的迭代步可能不够准确。

'fin-diff-grads'HessianMultiplyFcn 只有中等效率的原因是它们使用共轭梯度法。它们能准确估计目标函数的黑塞函数,但它们不生成最准确的迭代步。有关详细信息,请参阅fmincon 内点算法,及其关于求解公式 38 的 LDL 方法和共轭梯度法的讨论。

相关主题