使用符号数学与优化工具箱求解器

这个例子使用了:

打开实时脚本

这个例子展示了如何使用Symbolic Math Toolbox™函数雅可比矩阵而且matlabFunction为优化求解器提供解析导数。当您提供目标函数和约束函数的梯度和Hessians时，优化工具箱™求解器通常更准确和有效。

基于问题的优化可以自动计算和使用梯度;看到优化工具箱中的自动区分．有关使用自动区分的基于问题的示例，请参见基于优化变量的约束静电非线性优化．

在优化函数中使用符号计算有几个注意事项:

优化目标和约束函数应该用向量来定义x．然而，符号变量是标量值或复值，而不是向量值。这需要你在向量和标量之间进行转换。
优化梯度，有时是Hessians，应该在目标或约束函数的主体内计算。这意味着符号梯度或Hessian必须放置在目标或约束函数文件或函数句柄中的适当位置。
象征性地计算梯度和黑森可能很耗时。因此，您应该只执行此计算一次，并通过生成代码matlabFunction，在求解器执行期间调用。
方法计算符号表达式潜艇函数是耗时的。它的使用效率更高matlabFunction．
matlabFunction生成依赖于输入向量方向的代码。自fmincon使用列向量调用目标函数，调用时必须小心matlabFunction用符号变量的列向量。

第一个例子:Hessian的无约束最小化

最小化的目标函数为:

$f （ x_{1} ， x_{2} ）＝日志（ 1 + 3. {（ x_{2} - （ x_{1}^{3.} - x_{1} ））}^{2} + （ x_{1} - 4 / 3. ）^{2} ）．$

这个函数是正的，在处有一个唯一的最小值为零x1= 4/3,x2=(4/3)^3 - 4/3 = 1.0370…

我们把自变量写成x1而且x2因为在这种形式中，它们可以用作符号变量。作为向量的分量x它们会被写成x (1)而且x (2)．函数具有如下图所示的弯曲山谷。

信谊x1x2真正的X = [x1;x2];符号变量的列向量日志(f = 1 + 3 * (x2 - (x1 ^ 3 - x1)) ^ 2 + (x1 - 4/3) ^ 2)

f =
                  
                    $日志 （ {(x_{1} - \frac{4}{3.})}^{2} + 3. {(- {x_{1}}^{3.} + x_{1} + x_{2})}^{2} + 1 ）$

fsurf (f (2 - 2),“ShowContours”，“上”38)视图(127)

图中包含一个轴对象。axis对象包含一个functionsurface类型的对象。

计算f的梯度和Hessian:

Gradf =雅可比矩阵(f,x)。'% column gradf

gradf =
                  
                    $\begin{array}{l} （ \begin{array}{c} - \frac{6 (3. {x_{1}}^{2} - 1) (- {x_{1}}^{3.} + x_{1} + x_{2}) - 2 x_{1} + \frac{8}{3.}}{σ_{1}} \\ \frac{- 6 {x_{1}}^{3.} + 6 x_{1} + 6 x_{2}}{σ_{1}} \end{array} ） \\ 在哪里 \\ σ_{1} ＝ {(x_{1} - \frac{4}{3.})}^{2} + 3. {(- {x_{1}}^{3.} + x_{1} + x_{2})}^{2} + 1 \end{array}$

Hessf =雅可比矩阵(gradf,x)

hessf =
                  
                    $\begin{array}{l} （ \begin{array}{c} \frac{6 {(3. {x_{1}}^{2} - 1)}^{2} - 36 x_{1} (- {x_{1}}^{3.} + x_{1} + x_{2}) + 2}{σ_{2}} - \frac{{σ_{3.}}^{2}}{{σ_{2}}^{2}} & σ_{1} \\ σ_{1} & \frac{6}{σ_{2}} - \frac{{(- 6 {x_{1}}^{3.} + 6 x_{1} + 6 x_{2})}^{2}}{{σ_{2}}^{2}} \end{array} ） \\ 在哪里 \\ σ_{1} ＝ \frac{(- 6 {x_{1}}^{3.} + 6 x_{1} + 6 x_{2}) σ_{3.}}{{σ_{2}}^{2}} - \frac{18 {x_{1}}^{2} - 6}{σ_{2}} \\ σ_{2} ＝ {(x_{1} - \frac{4}{3.})}^{2} + 3. {(- {x_{1}}^{3.} + x_{1} + x_{2})}^{2} + 1 \\ σ_{3.} ＝ 6 (3. {x_{1}}^{2} - 1) (- {x_{1}}^{3.} + x_{1} + x_{2}) - 2 x_{1} + \frac{8}{3.} \end{array}$

的fminunc求解器期望传入一个向量x，并且SpecifyObjectiveGradient选项设置为真正的而且HessianFcn选项设置为“目标”，需要一个包含三个输出的列表:(f (x), gradf (x) hessf (x))．

matlabFunction从包含三个输入的列表中生成包含三个输出的列表。此外，使用var选项,matlabFunction接受矢量输入。

fh = matlabFunction(f,gradf,hessf，“var”, {x});

现在从点[-1,2]开始求解最小化问题:

选项= optimoptions(“fminunc”，.．.“SpecifyObjectiveGradient”,真的,.．.“HessianFcn”，“目标”，.．.“算法”，“信赖域”，.．.“显示”，“最后一次”）;[xfinal,fval,exitflag,output] = fminunc(fh，[-1;2]，options)

局部最小值。Fminunc停止，因为函数值相对于其初始值的最终变化小于函数公差值。

xfinal =2×11.3333 - 1.0370

Fval = 7.6623e-12

Exitflag = 3

输出=带字段的结构:迭代:14 funcCount: 15步长:0.0027 cgiterations: 11 firstorderopt: 3.4391e-05算法:'信任区域'消息:'本地最小可能....' constrviolation: []

将其与使用无梯度或Hessian信息的迭代次数进行比较。这需要“拟牛顿”算法。

选项= optimoptions(“fminunc”，“显示”，“最后一次”，“算法”，“拟牛顿”）;fh2 = matlabFunction(f，“var”, {x});% fh2 =目标没有梯度或黑森[xfinal,fval,exitflag,output2] = fminunc(fh2，[-1;2]，options)

找到局部极小值。优化完成，因为梯度的大小小于最优性公差的值。

xfinal =2×11.3333 - 1.0371

Fval = 2.1985e-11

Exitflag = 1

output2 =带字段的结构:迭代:18 funcCount: 81步长:2.1164e-04 lssteplth: 1 firstorderopt: 2.4587e-06算法:'准牛顿'消息:'本地最小值发现....'

当使用梯度和Hessians时，迭代次数会更少，并且函数计算也会大大减少:

sprintf (['有%d次迭代使用梯度'.．.“还有黑森语，但是没有。”]，.．.output.iterations output2.iterations)

ans =“使用梯度和Hessian的迭代有14次，但没有它们的迭代有18次。”

sprintf (['有%d个使用梯度的函数计算'.．.“还有黑森语，但是没有。”]，.．.output.funcCount output2.funcCount)

ans =“使用梯度和Hessian的函数计算有15个，但没有它们的有81个。”

第二个例子:使用fmincon内点算法的约束最小化

我们考虑相同的目标函数和起点，但现在有两个非线性约束:

$5 \sinh （ x_{2} / 5 ） \geq x_{1}^{4}$

$5 双曲正切（ x_{1} / 5 ） \geq x_{2}^{2} - 1 ．$

约束使优化远离全局最小值[1.333,1.037]。可视化这两个约束条件:

[X,Y] = meshgrid(-2:.01:3);Z = (5*sinh(y /5) >= x ^4);满足第一个约束时% Z=1，否则Z=0Z = Z+ 2*(5*tanh(X./5) >= Y.^2 - 1);% Z=2表示第二种满足，Z=3表示两种都满足冲浪(X, Y, Z,“线型”，“没有”）;FIG = gcf;fig.Color =' w '；%白色背景视图(0,90)在plot3(。43.96, .0373, 4,“o”，“MarkerEdgeColor”，“r”，“MarkerSize”8);最佳点包含(“x”) ylabel (“y”)举行从

图中包含一个轴对象。坐标轴对象包含两个类型为surface、line的对象。

我们在最佳点周围画了一个小红圈。

下面是可行区域上的目标函数图，该区域满足两个约束条件，如图中暗红色所示，在最优点周围有一个小红圈:

W = log(1 + 3*(Y - (X ^3 - X))。²+ (x - 4/3).²);% W =目标函数W(Z < 3) = nan;%仅在满足约束的地方绘图冲浪(X, Y, W,“线型”，“没有”）;视图(68年,20)在plot3(。43.96, .0373, .8152,“o”，“MarkerEdgeColor”，“r”，.．.“MarkerSize”8);最佳点包含(“x”) ylabel (“y”) zlabel (“z”)举行从

图中包含一个轴对象。坐标轴对象包含两个类型为surface、line的对象。

非线性约束必须写成这种形式C (x) <= 0．我们计算所有的符号约束及其导数，并将它们放在使用的函数句柄中matlabFunction．

约束的梯度应为列向量;它们必须以矩阵的形式放在目标函数中，矩阵的每一列表示一个约束函数的梯度。这是转置的形式由雅可比矩阵我们取下面的转置。

我们将非线性约束放入函数句柄中。fmincon期望非线性约束和梯度按顺序输出[c ceq gradc gradceq]．由于没有非线性等式约束，我们输出［］为量表信而且gradceq．

C1 = x1^4 - 5*sinh(x2/5);C2 = x2^2 - 5*tanh(x1/5) - 1;C = [c1 c2];Gradc =雅可比矩阵(c,x).';%转置到正确的形式约束= matlabFunction(c，[]，gradc，[]，“var”, {x});

内点算法要求其Hessian函数作为一个单独的函数来写，而不是作为目标函数的一部分。这是因为非线性约束函数需要在其黑森函数中包含这些约束。它的黑森是拉格朗日的黑森;更多信息请参见用户指南。

Hessian函数有两个输入参数:位置向量x，拉格朗日乘子结构lambda。结构中用于非线性约束的部分是lambda.ineqnonlin而且lambda.eqnonlin．对于当前的约束，没有线性等式，所以我们使用两个乘数lambda.ineqnonlin (1)而且lambda.ineqnonlin (2)．

在第一个例子中，我们计算了目标函数的Hessian。现在我们计算两个约束函数的Hessians，并使函数处理版本matlabFunction．

Hessc1 =雅可比矩阵(gradc(:，1)，x);% constraint =第一个c列Hessc2 =雅可比矩阵(gradc(:，2)，x);hessf = matlabFunction“var”, {x});hessc1h = matlabFunction(hessc1，“var”, {x});hessc2h = matlabFunction(hessc2，“var”, {x});

为了得到最终的黑森，我们将三个黑森放在一起，在约束函数中添加适当的拉格朗日乘子。

Myhess = @(x,lambda)(hessfh(x) +.．.lambda.ineqnonlin (1) * hessc1h (x) +.．.lambda.ineqnonlin (2) * hessc2h (x));

设置选项以使用内点算法、梯度和Hessian，让目标函数同时返回目标和梯度，并运行求解器:

选项= optimoptions(“fmincon”，.．.“算法”，“内点”，.．.“SpecifyObjectiveGradient”,真的,.．.“SpecifyConstraintGradient”,真的,.．.“HessianFcn”myhess,.．.“显示”，“最后一次”）;% fh2 =目标没有黑森fh2 = matlabFunction(f,gradf，“var”, {x});[xfinal,fval,exitflag,output] = fmincon(fh2，[-1;2]，.．.[],[],[],[],[],[], 约束,选项)

找到满足约束条件的局部最小值。优化完成是因为目标函数在可行方向上不递减，在最优性容差值范围内，约束条件满足在约束容差值范围内。

xfinal =2×10.4396 - 0.0373

Fval = 0.8152

Exitflag = 1

输出=带字段的结构:迭代:10 funcCount: 13 constrviolation: 0步长:1.9160e-06算法:'内部点' firstorderopt: 1.9217e-08 cgiterations: 0消息:'局部最小值发现满足约束....' bestviable: [1x1 struct]

同样，求解器在提供梯度和Hessian的情况下进行的迭代和函数计算比不提供梯度和Hessian时要少得多:

选项= optimoptions(“fmincon”，“算法”，“内点”，.．.“显示”，“最后一次”）;% fh3 =目标没有梯度或黑森fh3 = matlabFunction(f，“var”, {x});%无梯度约束:约束= matlabFunction(c，[]，“var”, {x});[xfinal,fval,exitflag,output2] = fmincon(fh3，[-1;2]，.．.[],[],[],[],[],[], 约束,选项)

找到了目标函数值较低的可行点。找到满足约束条件的局部最小值。优化完成是因为目标函数在可行方向上不递减，在最优性容差值范围内，约束条件满足在约束容差值范围内。

xfinal =2×10.4396 - 0.0373

Fval = 0.8152

Exitflag = 1

output2 =带字段的结构:迭代:18 funcCount: 57 constrviolation: 0步长:9.6632e-07算法:'内部点' firstorderopt: 3.8435e-07 cgiterations: 0消息:'局部最小值发现满足约束....' bestviable: [1x1 struct]

sprintf (['有%d次迭代使用梯度'.．.“还有黑森语，但是没有。”]，.．.output.iterations output2.iterations)

ans =“使用梯度和Hessian的迭代有10次，但没有它们的迭代有18次。”

sprintf (['有%d个使用梯度的函数计算'.．.“还有黑森语，但是没有。”]，.．.output.funcCount output2.funcCount)

ans =“使用梯度和Hessian的函数计算有13个，但没有它们的有57个。”

清理符号变量

本例中使用的符号变量被假定为实数。要从符号引擎工作区中清除这个假设，仅仅删除变量是不够的。您必须使用该语法清除变量的假设

假设((x1, x2),“清楚”）

当以下命令的输出为空时，所有的假设都将被清除

假设((x1, x2))

ans =空sym: 1乘0

使用符号数学与优化工具箱求解器

第一个例子:Hessian的无约束最小化

第二个例子:使用fmincon内点算法的约束最小化

清理符号变量

相关的话题