基于问题的工作流中的供应衍生工具——MATLAB和Simulink金宝appgydF4y2Ba

基于问题的工作流中的供应衍生品gydF4y2Ba

为什么包括衍生品吗?gydF4y2Ba

这个例子展示了如何在基于非线性问题的优化中包含导数信息，当自动导数不适用时，或者当你想要包含一个Hessian，它没有使用自动微分提供。在优化中包含梯度或Hessian可以带来以下好处:gydF4y2Ba

更健壮的结果。有限差分步骤有时会到达目标或非线性约束函数未定义、非有限或复杂的点。gydF4y2Ba
提高准确性。解析梯度比有限差分估计更精确。gydF4y2Ba
更快的收敛。加入Hessian可以导致更快的收敛，这意味着更少的迭代。gydF4y2Ba
提升的性能解析梯度的计算速度比有限差分估计要快，特别是对于稀疏结构的问题。然而，对于复杂的表达式，解析梯度的计算可能会比较慢。gydF4y2Ba

自动微分在最优化中的应用gydF4y2Ba

从R2020b开始gydF4y2Ba解决gydF4y2Ba函数可以对支持的函数使用自动微分，以便为求解器提供梯度。金宝app这些自动导数只适用于梯度(一阶导数)，而不适用于Hessians(二阶导数)。gydF4y2Ba

自动区分适用于您不使用gydF4y2Bafcn2optimexprgydF4y2Ba创建一个目标函数或约束函数。如果你需要的话gydF4y2Bafcn2optimexprgydF4y2Ba，这个例子展示了如何包含导数信息。gydF4y2Ba

在基于问题的优化中使用导数而不自动微分的方法是使用转换问题gydF4y2Baprob2structgydF4y2Ba，然后编辑生成的目标函数和约束函数。这个例子展示了一种混合方法，其中自动微分为部分目标函数提供导数。gydF4y2Ba

创建优化问题gydF4y2Ba

与控制变量gydF4y2BaxgydF4y2Ba和gydF4y2BaygydF4y2Ba，使用目标函数gydF4y2Ba

$\begin{array}{l} fun1gydF4y2Ba ＝gydF4y2Ba One hundred.gydF4y2Ba {（gydF4y2Ba ygydF4y2Ba -gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {（gydF4y2Ba 1gydF4y2Ba -gydF4y2Ba xgydF4y2Ba ）gydF4y2Ba}^{2gydF4y2Ba} \\ fun2gydF4y2Ba ＝gydF4y2Ba besseljgydF4y2Ba （gydF4y2Ba 1gydF4y2Ba ，gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba \\ 目标= fun1gydF4y2Ba +gydF4y2Ba fun2gydF4y2Ba ．gydF4y2Ba \end{array}$

包括的平方和的约束gydF4y2BaxgydF4y2Ba和gydF4y2BaygydF4y2Ba不超过4。gydF4y2Ba

fun2gydF4y2Ba不基于支持的函数进行优化表达式;金宝app看到gydF4y2Ba金宝app支持优化变量和表达式的操作gydF4y2Ba．因此,包括gydF4y2Bafun2gydF4y2Ba在优化问题中，必须使用以下方法将其转换为优化表达式gydF4y2Bafcn2optimexprgydF4y2Ba．gydF4y2Ba

要在支持的功能上使用AD，请设置不金宝app支持的功能的问题gydF4y2Bafun2gydF4y2Ba,包括gydF4y2Bafun2gydF4y2Ba以后。gydF4y2Ba

概率= optimproblem;x = optimvar (gydF4y2Ba“x”gydF4y2Ba）;y = optimvar (gydF4y2Ba“y”gydF4y2Ba）;Fun1 = 100*(y - x^2)^2 + (1 - x)^2;概率。目标= fun1;x^2 + y^2 <= 4;gydF4y2Ba

将问题转换为基于求解器的形式gydF4y2Ba

包括的导数gydF4y2Bafun2gydF4y2Ba，首先转换问题没有gydF4y2Bafun2gydF4y2Ba在结构中使用gydF4y2Baprob2structgydF4y2Ba．gydF4y2Ba

问题= prob2struct(概率,gydF4y2Ba．..gydF4y2Ba“ObjectiveFunctionName”gydF4y2Ba，gydF4y2Ba“generatedObjectiveBefore”gydF4y2Ba）;gydF4y2Ba

在转换期间,gydF4y2Baprob2structgydF4y2Ba创建表示目标和非线性约束函数的函数文件。默认情况下，这些函数都有名称gydF4y2Ba“generatedObjective.m”gydF4y2Ba和gydF4y2Ba“generatedConstraints.m”gydF4y2Ba．无目标函数文件gydF4y2Bafun2gydF4y2Ba是gydF4y2Ba“generatedObjectiveBefore.m”gydF4y2Ba．gydF4y2Ba

生成的目标函数和约束函数包括梯度。gydF4y2Ba

计算导数和跟踪变量gydF4y2Ba

计算相关的导数gydF4y2Bafun2gydF4y2Ba．如果您有Symbolic Math Toolbox™许可证，您可以使用gydF4y2Ba梯度gydF4y2Ba(符号数学工具箱)gydF4y2Ba或gydF4y2Ba雅可比矩阵gydF4y2Ba(符号数学工具箱)gydF4y2Ba函数来帮助计算导数。看到gydF4y2Ba使用符号数学工具箱™计算梯度和HessiansgydF4y2Ba．gydF4y2Ba

基于求解器的方法只有一个控制变量。每个优化变量(gydF4y2BaxgydF4y2Ba和gydF4y2BaygydF4y2Ba，在本例中)是控制变量的一部分。在生成的目标函数文件中可以找到优化变量与单个控制变量的映射关系gydF4y2Ba“generatedObjectiveBefore.m”gydF4y2Ba．文件的开头包含这些代码行或类似的行。gydF4y2Ba

% %变量指标。gydF4y2Baxidx = 1;yidx = 2;gydF4y2Ba将基于求解器的变量映射到基于问题的变量。gydF4y2Bax =数据源(xidx);y =数据源(yidx);gydF4y2Ba

在这段代码中，控制变量具有名称gydF4y2Ba数据源gydF4y2Ba．gydF4y2Ba

或者，您可以使用gydF4y2BavarindexgydF4y2Ba．gydF4y2Ba

idx = varindex(概率);disp (idx.x)gydF4y2Ba

1gydF4y2Ba

disp (idx.y)gydF4y2Ba

2gydF4y2Ba

整个目标函数包括gydF4y2Bafun2gydF4y2Ba：gydF4y2Ba

Fun2 = besselj(1,x^2 + y^2);gydF4y2Ba

使用标准微积分，计算gydF4y2Bagradfun2gydF4y2Ba，梯度gydF4y2Bafun2gydF4y2Ba．gydF4y2Ba

$gradfun2gydF4y2Ba ＝gydF4y2Ba ［gydF4y2Ba \begin{matrix} 2gydF4y2Ba xgydF4y2Ba （gydF4y2Ba besseljgydF4y2Ba （gydF4y2Ba 0gydF4y2Ba ，gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba -gydF4y2Ba besseljgydF4y2Ba （gydF4y2Ba 1gydF4y2Ba ，gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba /gydF4y2Ba （gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba ）gydF4y2Ba \\ 2gydF4y2Ba ygydF4y2Ba （gydF4y2Ba besseljgydF4y2Ba （gydF4y2Ba 0gydF4y2Ba ，gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba -gydF4y2Ba besseljgydF4y2Ba （gydF4y2Ba 1gydF4y2Ba ，gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba /gydF4y2Ba （gydF4y2Ba {xgydF4y2Ba}^{2gydF4y2Ba} +gydF4y2Ba {ygydF4y2Ba}^{2gydF4y2Ba} ）gydF4y2Ba ）gydF4y2Ba \end{matrix} ］gydF4y2Ba ．gydF4y2Ba$

编辑目标和约束文件gydF4y2Ba

编辑gydF4y2Ba“generatedObjectiveBefore.m”gydF4y2Ba包括gydF4y2Bafun2gydF4y2Ba．gydF4y2Ba

%%计算目标函数。Arg1 = (y - x.^2);最长= 100;长度= __arg1。^ 2;Arg4 = (1 - x);Obj = ((arg2 .* arg3) + arg4.^2);gydF4y2BaSSQ = x^2 + y^2;fun2 = besselj (ssq);Obj = Obj + fun2;gydF4y2Ba

通过编辑在目标函数文件中包含计算出的梯度gydF4y2Ba“generatedObjectiveBefore.m”gydF4y2Ba．如果你有一个软件版本，不执行梯度计算，包括所有这些线。如果您的软件执行梯度计算，只包括粗体线，它计算的梯度gydF4y2Bafun2gydF4y2Ba．gydF4y2Ba

计算客观梯度。If nargout > 1 arg5 = 1;Arg6 = 0 ([2, 1]);arg6 (xidx:) = (- (arg5。* 2 * (arg4 (:)))) + ((-(( arg5。*最长(:))。* 2 * (__arg1(:))))。* 2 * (x (:)));arg6 (yidx:) = ((arg5。*最长(:))。* 2。* (__arg1 (:)));研究生= arg6 (:);gydF4y2Baarg7 = besselj(0、ssq);Arg8 = arg7 - fun2/ssq;gfun = (2 * x * arg8;…2 * y * arg8];Grad = Grad + gfun;gydF4y2Ba结束gydF4y2Ba

回想一下非线性约束是gydF4y2BaX ^2 + y^2 <= 4gydF4y2Ba．这个约束函数的梯度是gydF4y2Ba2 * (x, y)gydF4y2Ba．如果您的软件计算约束梯度并将其包含在生成的约束文件中，那么您不需要再做任何事情。如果你的软件不计算约束梯度，那么通过编辑包括非线性约束的梯度gydF4y2Ba“generatedConstraints.m”gydF4y2Ba包括这些行。gydF4y2Ba

%%此处插入渐变计算。如果你包含一个梯度，通过设置% specificconstraintgradient选项为true通知求解器。如果nargout > 2gydF4y2BacineqGrad = 2 * (x, y);gydF4y2BaceqGrad = [];结束gydF4y2Ba

运行问题有和没有梯度gydF4y2Ba

同时使用基于求解器和基于问题(无梯度)的方法运行问题，以查看差异。要使用衍生信息运行基于求解器的问题，请在问题结构中创建适当的选项。gydF4y2Ba

选择= optimoptions (gydF4y2Ba“fmincon”gydF4y2Ba，gydF4y2Ba“SpecifyObjectiveGradient”gydF4y2Ba,真的,gydF4y2Ba．..gydF4y2Ba“SpecifyConstraintGradient”gydF4y2Ba,真正的);问题。选择=选项;gydF4y2Ba

非线性问题需要一个非空gydF4y2Bax0gydF4y2Ba字段在问题结构中。gydF4y2Ba

x0 = (1, 1);问题。x0 = x0;gydF4y2Ba

调用gydF4y2BafmincongydF4y2Ba关于问题结构。gydF4y2Ba

[xsolver, fvalsolver exitflagsolver outputsolver] = fmincon(问题)gydF4y2Ba

找到满足约束条件的局部最小值。优化完成是因为目标函数在可行方向上不减小到最优性公差的值内，约束条件满足到约束公差的值内。<停止条件详细信息> xsolver = 1.2494 1.5617 fvalsolver = -0.0038 exitflagsolver = 1 outputsolver = struct with fields: iterations: 15 funcCount: 32 construct: 0 stepsize: 1.5569e-06 algorithm: 'interior-point' firstorderopt: 2.2058e-08 cgiterations: 7 message: '↵Local minimum found that满足约束。↵↵优化完成，因为目标函数在↵可行方向上不减小到最优性容差值内，↵和约束条件满足到约束容差值内。↵↵<停止条件详细信息>↵↵优化完成:相对一阶优化测度，2.125635e-08，↵小于选项。OptimalityTolerance = 1.000000e-06，相对最大约束↵违例，0.000000e+00，小于选项。ConstraintTolerance = 1.000000 e-06。↵↵'最佳可行:[1×1 struct]gydF4y2Ba

将这个解与由gydF4y2Ba解决gydF4y2Ba没有导数信息。gydF4y2Ba

init。x＝x0（1）;init。y＝x0（2）;F2 = @(x,y)besselj(1,x^2 + y^2);fun2 = fcn2optimexpr (f2, x, y); prob.Objective = prob.Objective + fun2; [xproblem,fvalproblem,exitflagproblem,outputproblem] = solve(prob,init,．..gydF4y2Ba“ConstraintDerivative”gydF4y2Ba，gydF4y2Ba“有限差分”gydF4y2Ba，gydF4y2Ba．..gydF4y2Ba“ObjectiveDerivative”gydF4y2Ba，gydF4y2Ba“有限差分”gydF4y2Ba）gydF4y2Ba

使用fmincon解决问题。找到满足约束条件的局部最小值。优化完成是因为目标函数在可行方向上不减小到最优性公差的值内，约束条件满足到约束公差的值内。 xproblem = struct with fields: x: 1.2494 y: 1.5617 fvalproblem = -0.0038 exitflagproblem = OptimalSolution outputproblem = struct with fields: iterations: 15 funcCount: 64 construction: 0 stepsize: 1.5571e-06 algorithm: ' internal -point' firstderopt: 6.0139e-07 cgiterations: 7 message:↵局部最小值满足约束条件。↵↵优化完成，因为目标函数在↵可行方向上不减小到最优性容差值内，↵和约束条件满足到约束容差值内。↵↵<停止条件详细信息>↵↵优化完成:相对一阶优化测度，5.795368e-07，↵小于选项。OptimalityTolerance = 1.000000e-06，相对最大约束↵违例，0.000000e+00，小于选项。ConstraintTolerance = 1.000000 e-06。↵↵'最佳可行:[1×1 struct]objectivederivative: "finite-differences" constraintderivative: "closed-form" solver: 'fmincon'

这两个解决金宝搏官方网站方案报告相同的函数值以显示精度，并且需要相同的迭代次数。然而，有梯度的解决方案只需要32个函数计算，而没有梯度的解决方案需要64个函数计算。gydF4y2Ba

包括黑森gydF4y2Ba

要包括一个黑森，你必须使用gydF4y2Baprob2structgydF4y2Ba，即使您的所有函数都支持优化表达式。金宝app这个例子展示了如何使用HessiangydF4y2BafmincongydF4y2Ba内点gydF4y2Ba算法。的gydF4y2BafminuncgydF4y2Ba信赖域gydF4y2Ba算法和gydF4y2BafmincongydF4y2Batrust-region-reflectivegydF4y2Ba算法使用不同的语法;看到gydF4y2BaHessian用于fminunc信任区域或fmincon信任区域反射算法gydF4y2Ba．gydF4y2Ba

中描述的gydF4y2Bafmincon内点算法的HessiangydF4y2Ba，黑森就是拉格朗日的黑森。gydF4y2Ba

{\nablagydF4y2Ba}_{xgydF4y2Ba xgydF4y2Ba}^{2gydF4y2Ba} lgydF4y2Ba （gydF4y2Ba xgydF4y2Ba ，gydF4y2Ba λgydF4y2Ba ）gydF4y2Ba ＝gydF4y2Ba {\nablagydF4y2Ba}^{2gydF4y2Ba} fgydF4y2Ba （gydF4y2Ba xgydF4y2Ba ）gydF4y2Ba +gydF4y2Ba \sumgydF4y2Ba {λgydF4y2Ba}_{ggydF4y2Ba ，gydF4y2Ba 我gydF4y2Ba} {\nablagydF4y2Ba}^{2gydF4y2Ba} {ggydF4y2Ba}_{我gydF4y2Ba} （gydF4y2Ba xgydF4y2Ba ）gydF4y2Ba +gydF4y2Ba \sumgydF4y2Ba {λgydF4y2Ba}_{hgydF4y2Ba ，gydF4y2Ba 我gydF4y2Ba} {\nablagydF4y2Ba}^{2gydF4y2Ba} {hgydF4y2Ba}_{我gydF4y2Ba} （gydF4y2Ba xgydF4y2Ba ）gydF4y2Ba ．gydF4y2Ba

（1）gydF4y2Ba

通过创建函数文件来计算Hessian，并创建gydF4y2BaHessianFcngydF4y2Ba选择gydF4y2BafmincongydF4y2Ba使用黑森枪。为了在这种情况下创建Hessian，分别创建目标和非线性约束的二阶导数。gydF4y2Ba

本例中目标函数的二阶导数矩阵比较复杂。其目标函数列表gydF4y2Bahessianfun (x)gydF4y2Ba是由符号数学工具箱创建的，使用的方法与gydF4y2Ba使用符号数学工具箱™计算梯度和HessiansgydF4y2Ba．gydF4y2Ba

函数gydF4y2Ba高频= hessfun (in)gydF4y2Ba% HESSFUNgydF4y2Ba% hf = hessfun (in1)gydF4y2Ba这个函数是由symbol Math Toolbox 8.6版本生成的。gydF4y2Ba% 10 - 8月- 2020 10:50:44gydF4y2Bax = in (1:);y = in (2:);t2 = x ^ 2;t4 = y ^ 2;t6 = x。* 4.0 e + 2;t3 = t2。^ 2;t5 = t4。^ 2;t7 = t4;t8 = t6;t9 = t2 + t4; t10 = t2.*t4.*2.0; t11 = besselj(0,t9); t12 = besselj(1,t9); t13 = t2+t7; t14 = 1.0./t9; t16 = t3+t5+t10-2.0; t15 = t14.^2; t17 = t11.*t14.*x.*y.*4.0; t19 = t11.*t13.*t14.*2.0; t18 = -t17; t20 = t12.*t15.*t16.*x.*y.*4.0; t21 = -t20; t22 = t8+t18+t21; hf = reshape([t2.*1.2e+3-t19-y.*4.0e+2-t12.*t15.*．..gydF4y2Ba(t2 * -3.0 + t4 + t2。* t5。* 2.0 + t3, t4 *。* 4.0 + t2。^ 3。* 2.0)* 2.0 + 2.0,gydF4y2Ba．..gydF4y2Bat22 t22,gydF4y2Ba．..gydF4y2Bat19-t12。* * (t2-t4 t15。* 3.0 + t2 * t5。* 4.0 +gydF4y2Ba．..gydF4y2Bat3, t4 *。* 2.0 + t4。^ 3 * 2.0)* 2.0 + 2.0 e + 2], (2, 2));gydF4y2Ba

相反，非线性不等式约束的Hessian是简单的;它是单位矩阵的两倍。gydF4y2Ba

hessianc = 2 *眼(2);gydF4y2Ba

创建作为函数柄的拉格朗日的Hessian。gydF4y2Ba

H = @(x,lam)(hessianfun(x) + hessianc*lam.ineqnonlin);gydF4y2Ba

创建使用此Hessian的选项。gydF4y2Ba

problem.options.HessianFcn = H;gydF4y2Ba

解决问题并显示迭代次数和函数计算次数。解和之前差不多。gydF4y2Ba

[xhess, fvalhess exitflaghess outputhess] = fmincon(问题);disp (outputhess.iterations) disp (outputhess.funcCount)gydF4y2Ba

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

这一次,gydF4y2BafmincongydF4y2Ba只需要8次迭代，而不是15次，只需要10次函数计算，而不是32次。总之，提供一个解析Hessian计算可以提高求解过程的效率，但发展一个函数来计算Hessian可能是困难的。gydF4y2Ba

另请参阅gydF4y2Ba

fcn2optimexprgydF4y2Ba|gydF4y2Baprob2structgydF4y2Ba|gydF4y2BavarindexgydF4y2Ba