## Comparison of Six Solvers

### Function to Optimize

This example shows how to minimize Rastrigin’s function with six solvers. Each solver has its own characteristics. The characteristics lead to different solutions and run times. The results, examined in Compare Syntax and Solutions, can help you choose an appropriate solver for your own problems.

Rastrigin’s function has many local minima, with a global minimum at (0,0):

`$\text{Ras}\left(x\right)=20+{x}_{1}^{2}+{x}_{2}^{2}-10\left(\mathrm{cos}2\pi {x}_{1}+\mathrm{cos}2\pi {x}_{2}\right).$`

Usually you don't know the location of the global minimum of your objective function. To show how the solvers look for a global solution, this example starts all the solvers around the point `[20,30]`, which is far from the global minimum.

The `rastriginsfcn.m` file implements Rastrigin’s function. This file comes with Global Optimization Toolbox software. This example employs a scaled version of Rastrigin’s function with larger basins of attraction. For information, see Basins of Attraction.

`rf2 = @(x)rastriginsfcn(x/10);` This example minimizes `rf2` using the default settings of `fminunc` (an Optimization Toolbox™ solver), `patternsearch`, and `GlobalSearch`. The example also uses `ga` and `particleswarm` with nondefault options to start with an initial population around the point `[20,30]`. Because `surrogateopt` requires finite bounds, the example uses `surrogateopt` with lower bounds of `-70` and upper bounds of `130` in each variable.

### Six Solution Methods

#### fminunc

To solve the optimization problem using the `fminunc` Optimization Toolbox solver, enter:

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

`fminunc` returns

```Local minimum found. Optimization completed because the size of the gradient is less than the default value of the function tolerance. xf = 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` is the minimizing point.

• `ff` is the value of the objective, `rf2`, at `xf`.

• `flf` is the exit flag. An exit flag of `1` indicates `xf` is a local minimum.

• `of` is the output structure, which describes the `fminunc` calculations leading to the solution.

#### patternsearch

To solve the optimization problem using the `patternsearch` Global Optimization Toolbox solver, enter:

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

`patternsearch` returns

```Optimization terminated: mesh size less than options.MeshTolerance. xp = 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: 'Optimization terminated: mesh size less than options.MeshTolerance.'```
• `xp` is the minimizing point.

• `fp` is the value of the objective, `rf2`, at `xp`.

• `flp` is the exit flag. An exit flag of `1` indicates `xp` is a local minimum.

• `op` is the output structure, which describes the `patternsearch` calculations leading to the solution.

#### ga

To solve the optimization problem using the `ga` Global Optimization Toolbox solver, enter:

```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); [xga,fga,flga,oga] = ga(rf2,2,[],[],[],[],[],[],[],opts)```

`initpop` is a 20-by-2 matrix. Each row of `initpop` has mean `[20,30]`, and each element is normally distributed with standard deviation `10`. The rows of `initpop` form an initial population matrix for the `ga` solver.

`opts` is the options that set `initpop` as the initial population.

The final line calls `ga`, using the options.

`ga` uses random numbers, and produces a random result. In this case `ga` returns:

```Optimization terminated: maximum number of generations exceeded. xga = -0.0042 -0.0024 fga = 4.7054e-05 flga = 0 oga = struct with fields: problemtype: 'unconstrained' rngstate: [1×1 struct] generations: 200 funccount: 9453 message: 'Optimization terminated: maximum number of generations exceeded.' maxconstraint: []```
• `xga` is the minimizing point.

• `fga` is the value of the objective, `rf2`, at `xga`.

• `flga` is the exit flag. An exit flag of `0` indicates that `ga` reached a function evaluation limit or an iteration limit. In this case, `ga` reached an iteration limit.

• `oga` is the output structure, which describes the `ga` calculations leading to the solution.

#### particleswarm

Like `ga`, `particleswarm` is a population-based algorithm. So for a fair comparison of solvers, initialize the particle swarm to the same population as `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 = 9.9496 0.0000 fpso = 0.9950 flgpso = 1 opso = struct with fields: rngstate: [1×1 struct] iterations: 56 funccount: 1140 message: 'Optimization ended: relative change in the objective value ↵over the last OPTIONS.MaxStallIterations iterations is less than OPTIONS.FunctionTolerance.'```
• `xpso` is the minimizing point.

• `fpso` is the value of the objective, `rf2`, at `xpso`.

• `flgpso` is the exit flag. An exit flag of `1` indicates `xpso` is a local minimum.

• `opso` is the output structure, which describes the `particleswarm` calculations leading to the solution.

#### surrogateopt

`surrogateopt` does not require a start point, but does require finite bounds. Set bounds of –70 to 130 in each component. To have the same sort of output as the other solvers, disable the default plot function.

```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 = -0.0033 0.0005 fsur = 2.2456e-05 flgsur = 0 osur = struct with fields: elapsedtime: 2.3877 funccount: 200 rngstate: [1×1 struct] message: 'Surrogateopt stopped because it exceeded the function evaluation limit set by ↵'options.MaxFunctionEvaluations'.' ```
• `xsur` is the minimizing point.

• `fsur` is the value of the objective, `rf2`, at `xsur`.

• `flgsur` is the exit flag. An exit flag of `0` indicates that `surrogateopt` halted because it ran out of function evaluations or time.

• `osur` is the output structure, which describes the `surrogateopt` calculations leading to the solution.

#### GlobalSearch

To solve the optimization problem using the `GlobalSearch` solver, enter:

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

`problem` is an optimization problem structure. `problem` specifies the `fmincon` solver, the `rf2` objective function, and `x0=[20,30]`. For more information on using `createOptimProblem`, see Create Problem Structure.

Note

You must specify `fmincon` as the solver for `GlobalSearch`, even for unconstrained problems.

`gs` is a default `GlobalSearch` object. The object contains options for solving the problem. Calling `run(gs,problem)` runs `problem` from multiple start points. The start points are random, so the following result is also random.

In this case, the run returns:

```GlobalSearch stopped because it analyzed all the trial points. All 10 local solver runs converged with a positive local solver exit flag. xg = 1.0e-07 * -0.1405 -0.1405 fg = 0 flg = 1 og = struct with fields: funcCount: 2350 localSolverTotal: 10 localSolverSuccess: 10 localSolverIncomplete: 0 localSolverNoSolution: 0 message: 'GlobalSearch stopped because it analyzed all the trial points.↵↵All 10 local solver runs converged with a positive local solver exit flag.' ```
• `xg` is the minimizing point.

• `fg` is the value of the objective, `rf2`, at `xg`.

• `flg` is the exit flag. An exit flag of `1` indicates all `fmincon` runs converged properly.

• `og` is the output structure, which describes the `GlobalSearch` calculations leading to the solution.

### Compare Syntax and Solutions

One solution is better than another if its objective function value is smaller than the other. The following table summarizes the results, accurate to one decimal.

ResultsfminuncpatternsearchgaparticleswarmsurrogateoptGlobalSearch
solution`[19.9 29.9]``[19.9 -9.9]``[0 0]``[10 0]``[0 0]``[0 0]`
objective`12.9``5``0``1``0``0`
# Fevals`15``174``9453``1140``200``2178`

These results are typical:

• `fminunc` quickly reaches the local solution within its starting basin, but does not explore outside this basin at all. `fminunc` has a simple calling syntax.

• `patternsearch` takes more function evaluations than `fminunc`, and searches through several basins, arriving at a better solution than `fminunc`. The `patternsearch` calling syntax is the same as that of `fminunc`.

• `ga` takes many more function evaluations than `patternsearch`. By chance it arrived at a better solution. In this case, `ga` found a point near the global optimum. `ga` is stochastic, so its results change with every run. `ga` has a simple calling syntax, but there are extra steps to have an initial population near `[20,30]`.

• `particleswarm` takes fewer function evaluations than `ga`, but more than `patternsearch`. In this case, `particleswarm` found a point with lower objective function value than `patternsearch`, but higher than `ga`. Because `particleswarm` is stochastic, its results change with every run. `particleswarm` has a simple calling syntax, but there are extra steps to have an initial population near `[20,30]`.

• `surrogateopt` stops when it reaches a function evaluation limit, which by default is 200 for a two-variable problem. `surrogateopt` has a simple calling syntax, but requires finite bounds. `surrogateopt` attempts to find a global solution, and in this case succeeded. Each function evaluation in `surrogateopt` takes a longer time than in most other solvers, because `surrogateopt` performs many auxiliary computations as part of its algorithm.

• `GlobalSearch` `run` takes the same order of magnitude of function evaluations as `ga` and `particleswarm`, searches many basins, and arrives at a good solution. In this case, `GlobalSearch` found the global optimum. Setting up `GlobalSearch` is more involved than setting up the other solvers. As the example shows, before calling `GlobalSearch`, you must create both a `GlobalSearch` object (`gs` in the example), and a problem structure (`problem`). Then, you call the `run` method with `gs` and `problem`. For more details on how to run `GlobalSearch`, see Workflow for GlobalSearch and MultiStart.

Watch now