六种求解器的比较

要优化的函数

此示例说明如何使用六个求解器来最小化拉斯特里金函数。每个求解器都有自己的特点。这些特点导致不同的解和运行时间。在中检查的结果可以帮助您为自己的问题选择合适的求解器。

拉斯特里金函数具有许多局部极小值，全局最小值位于 (0,0)。该函数定义为 $R a s (x)$ ：

$R a s (x) = 20 + x_{1}^{2} + x_{2}^{2} - 10 (\cos 2 π x_{1} + \cos 2 π x_{2}) .$

运行此示例时，可以使用计算拉斯特里金函数值的 rastriginsfcn.m 文件。此示例采用了拉斯特里金函数的缩放版本，具有更大的吸引盆。有关信息，请参阅吸引盆。创建缩放函数的曲面图。

rf2 = @(x)rastriginsfcn(x/10);
rf3 = @(x,y)reshape(rastriginsfcn([x(:)/10,y(:)/10]),size(x));
fsurf(rf3,[-30 30],'ShowContours','on')
title('rastriginsfcn([x/10,y/10])')
xlabel('x')
ylabel('y')

通常，您不知道目标函数的全局最小值的位置。为了展示求解器如何寻找全局解，此示例从点 [20,30] 附近启动所有求解器，该点远离全局最小值。

此示例使用 fminunc（Optimization Toolbox™ 求解器）、patternsearch 和 GlobalSearch 的默认设置最小化 rf2。该示例还使用具有非默认选项的 ga 和 particleswarm 来从点 [20,30] 周围的初始种群开始。因为 surrogateopt 需要有限边界，所以该示例使用 surrogateopt，每个变量的下界为 -70，上界为 130。

六种解决方法

`fminunc`

要使用 fminunc Optimization Toolbox 求解器解决优化问题，请输入：

rf2 = @(x)rastriginsfcn(x/10); % objective
x0 = [20,30]; % start point away from the minimum
[xf,ff,flf,of] = fminunc(rf2,x0)

Local minimum found.

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

xf = 1×2

   19.8991   29.8486

ff = 12.9344

flf = 1

of = struct with fields:
       iterations: 3
        funcCount: 15
         stepsize: 1.7776e-06
     lssteplength: 1
    firstorderopt: 5.9907e-09
        algorithm: 'quasi-newton'
          message: 'Local minimum found....'

xf 是最小化点。
ff 是 xf 处目标 rf2 的值。
flf 是退出标志。退出标志为 1 表示 xf 是局部最小值。
of 是输出结构体，它描述了导致解的 fminunc 计算。

`patternsearch`

要使用 patternsearch Global Optimization Toolbox 求解器解决优化问题，请输入：

rf2 = @(x)rastriginsfcn(x/10); % objective
x0 = [20,30]; % start point away from the minimum
[xp,fp,flp,op] = patternsearch(rf2,x0)

patternsearch stopped because the mesh size was less than options.MeshTolerance.

xp = 1×2

   19.8991   -9.9496

fp = 4.9748

flp = 1

op = struct with fields:
         function: @(x)rastriginsfcn(x/10)
      problemtype: 'unconstrained'
       pollmethod: 'gpspositivebasis2n'
    maxconstraint: []
     searchmethod: []
       iterations: 48
        funccount: 174
         meshsize: 9.5367e-07
         rngstate: [1x1 struct]
          message: 'patternsearch stopped because the mesh size was less than options.MeshTolerance.'

xp 是最小化点。
fp 是 xp 处目标 rf2 的值。
flp 是退出标志。退出标志为 1 表示 xp 是局部最小值。
op 是输出结构体，它描述了导致解的 patternsearch 计算。

`ga`

要使用 ga Global Optimization Toolbox 求解器解决优化问题，请输入：

rng default % for reproducibility
rf2 = @(x)rastriginsfcn(x/10); % objective
x0 = [20,30]; % start point away from the minimum
initpop = 10*randn(20,2) + repmat(x0,20,1);
opts = optimoptions('ga','InitialPopulationMatrix',initpop);

initpop 是一个 20×2 的矩阵。initpop 的每一行都有平均值 [20,30]，并且每个元素都呈标准差为 10 的正态分布。initpop 的行形成 ga 求解器的初始种群矩阵。
opts 是将 initpop 设置为初始种群的选项。
ga 使用随机数，并产生随机结果。
下一行使用选项调用 ga。

[xga,fga,flga,oga] = ga(rf2,2,[],[],[],[],[],[],[],opts)

ga stopped because it exceeded options.MaxGenerations.

xga = 1×2

   -0.0042   -0.0024

fga = 4.7054e-05

flga = 0

oga = struct with fields:
      problemtype: 'unconstrained'
         rngstate: [1x1 struct]
      generations: 200
        funccount: 9453
          message: 'ga stopped because it exceeded options.MaxGenerations.'
    maxconstraint: []
       hybridflag: []

xga 是最小化点。
fga 是 xga 处目标 rf2 的值。
flga 是退出标志。退出标志为 0 表示 ga 达到函数评估限制或迭代限制。在这种情况下，ga 达到了迭代极限。
oga 是输出结构体，它描述了导致解的 ga 计算。

`particleswarm`

与 ga 一样，particleswarm 是一种基于种群的算法。因此，为了公平地比较求解器，将粒子群初始化为与 ga 相同的种群。

rng default % for reproducibility
rf2 = @(x)rastriginsfcn(x/10); % objective
opts = optimoptions('particleswarm','InitialSwarmMatrix',initpop);
[xpso,fpso,flgpso,opso] = particleswarm(rf2,2,[],[],opts)

Optimization ended: relative change in the objective value 
over the last OPTIONS.MaxStallIterations iterations is less than OPTIONS.FunctionTolerance.

xpso = 1×2

    9.9496    0.0000

fpso = 0.9950

flgpso = 1

opso = struct with fields:
      rngstate: [1x1 struct]
    iterations: 56
     funccount: 1140
       message: 'Optimization ended: relative change in the objective value ...'
    hybridflag: []

`surrogateopt`

surrogateopt 不需要起点，但需要有限的边界。在每个分量中设置-70 到 130 的边界。为了获得与其他求解器相同类型的输出，请禁用默认绘图函数。

rng default % for reproducibility
lb = [-70,-70];
ub = [130,130];
rf2 = @(x)rastriginsfcn(x/10); % objective
opts = optimoptions('surrogateopt','PlotFcn',[]);
[xsur,fsur,flgsur,osur] = surrogateopt(rf2,lb,ub,opts)

surrogateopt stopped because it exceeded the function evaluation limit set by 
'options.MaxFunctionEvaluations'.

xsur = 1×2

    9.9494   -9.9502

fsur = 1.9899

flgsur = 0

osur = struct with fields:
        elapsedtime: 3.4809
          funccount: 200
    constrviolation: 0
               ineq: [1x0 double]
           rngstate: [1x1 struct]
            message: 'surrogateopt stopped because it exceeded the function evaluation limit set by ...'

xsur 是最小化点。
fsur 是 xsur 处目标 rf2 的值。
flgsur 是退出标志。退出标志为 0 表示 surrogateopt 因为函数计算或时间耗尽而停止。
osur 是输出结构体，它描述了导致解的 surrogateopt 计算。

`GlobalSearch`

要使用 GlobalSearch 求解器解决优化问题，请输入：

rf2 = @(x)rastriginsfcn(x/10); % objective
x0 = [20,30]; % start point away from the minimum
problem = createOptimProblem('fmincon','objective',rf2,...
    'x0',x0);
gs = GlobalSearch;

problem 是一个优化问题结构体。problem 指定 fmincon 求解器、rf2 目标函数和 x0=[20,30]。有关使用 createOptimProblem 的详细信息，请参阅创建问题结构体。

注意：即使对于无约束的问题，您也必须指定 fmincon 作为 GlobalSearch 的求解器。

gs 是默认的 GlobalSearch 对象。该对象包含解决问题的选项。调用 run(gs,problem) 从多个起点运行问题。起点是随机的，所以接下来的结果也是随机的。

[xg,fg,flg,og] = run(gs,problem)

GlobalSearch stopped because it analyzed all the trial points.

All 7 local solver runs converged with a positive local solver exit flag.

xg = 1×2
10^-7 ×

   -0.1405   -0.1405

fg = 0

flg = 1

og = struct with fields:
                funcCount: 2217
         localSolverTotal: 7
       localSolverSuccess: 7
    localSolverIncomplete: 0
    localSolverNoSolution: 0
                  message: 'GlobalSearch stopped because it analyzed all the trial points....'

xg 是最小化点。
fg 是 xg 处目标 rf2 的值。
flg 是退出标志。退出标志为 1 表示所有 fmincon 运行均正确收敛。
og 是输出结构体，它描述了导致解的 GlobalSearch 计算。

比较语法和解决方案

如果一个解决方案的目标函数值小于另一个解决方案，则该解比另一个解决方案更好。下表总结了结果。

sols = [xf;
    xp;
    xga;
    xpso;
    xsur;
    xg];
fvals = [ff;
    fp;
    fga;
    fpso;
    fsur;
    fg];
fevals = [of.funcCount;
    op.funccount;
    oga.funccount;
    opso.funccount;
    osur.funccount;
    og.funcCount];
stats = table(sols,fvals,fevals);
stats.Properties.RowNames = ["fminunc" "patternsearch" "ga" "particleswarm" "surrogateopt" "GlobalSearch"];
stats.Properties.VariableNames = ["Solution" "Objective" "# Fevals"];
disp(stats)

                              Solution             Objective     # Fevals
                     __________________________    __________    ________

    fminunc               19.899         29.849        12.934        15  
    patternsearch         19.899        -9.9496        4.9748       174  
    ga                -0.0042178     -0.0024347    4.7054e-05      9453  
    particleswarm         9.9496       6.75e-07       0.99496      1140  
    surrogateopt          9.9494        -9.9502        1.9899       200  
    GlobalSearch     -1.4046e-08    -1.4046e-08             0      2217

这些结果是典型的：

fminunc 很快到达其起始盆地内的局部解，但根本没有探索该盆地之外的区域。fminunc 具有简单的调用语法。
patternsearch 比 fminunc 进行了更多的函数计算，并通过多个盆地进行了搜索，最终得出了比 fminunc 更好的解。patternsearch 调用语法与 fminunc 相同。
ga 比 patternsearch 需要更多的函数计算。它偶然得出了一个更好的解。在这种情况下，ga 找到一个接近全局最优的点。ga 是随机的，因此其结果会随着每次运行而变化。ga 有一个简单的调用语法，但需要额外的步骤来在 [20,30] 附近有一个初始种群。
particleswarm 所需的函数计算比 ga 少，但比 patternsearch 多。在这种情况下，particleswarm 找到一个目标函数值低于 patternsearch 但高于 ga 的点。因为 particleswarm 是随机的，所以它的结果每次运行都会发生变化。particleswarm 有一个简单的调用语法，但需要额外的步骤才能在 [20,30] 附近有一个初始种群。
surrogateopt 在达到函数评估限制时停止，对于双变量问题，默认值为 200。surrogateopt 具有简单的调用语法，但需要有限边界。surrogateopt 尝试寻找全局解，但在这种情况下没有成功。surrogateopt 中的每个函数评估比大多数其他求解器花费的时间更长，因为 surrogateopt 作为其算法的一部分执行许多辅助计算。
GlobalSearch run 与 ga 和 particleswarm 采用相同数量级的函数计算值，搜索了许多盆地，并找到了一个好的解。在这种情况下，GlobalSearch 找到了全局最优。设置 GlobalSearch 比设置其他求解器更复杂。如例所示，在调用 GlobalSearch 之前，必须创建一个 GlobalSearch 对象（示例中为 gs）和一个问题结构体（problem）。然后，使用 gs 和 problem 调用 run 方法。有关如何运行 GlobalSearch 的更多详细信息，请参阅 GlobalSearch 和 MultiStart 的工作流。