Search for Global Minimum Using patternsearch

Open Live Script

The patternsearch solver does not always converge to a global minimum. To increase your chances of obtaining a global minimum, restart patternsearch multiple times from random start points.

Problem with Multiple Local Minima

The sawtoothxy function has multiple local minima, and a global minimum at [0,0] with function value 0.

function f = sawtoothxy(x,y)
[t,r] = cart2pol(x,y); % change to polar coordinates
h = cos(2*t - 1/2)/2 + cos(t) + 2;
g = (sin(r) - sin(2*r)/2 + sin(3*r)/3 - sin(4*r)/4 + 4) ...
    .*r.^2./(r+1);
f = g.*h;
end

Create asymmetric linear inequality constraints for the problem, and a feasible initial point far from the global minimum.

A = [1 1
    1 -1
    -1 1
    -1 -1];
b = [200,200,20,10];
fun = @(t)sawtoothxy(t(1),t(2));
rng(10) % Sets random number stream
x0 = [100,-50] + 10*randn(1,2);

Check that the initial point is feasible, and evaluate the sawtoothxy function at that point.

A*x0' - b' % Should be negative for a feasible initial point.
ans = 4×1

 -158.1931
  -28.8906
 -191.1094
  -51.8069


fun(x0)
ans = 
1.1681e+03

The patternsearch "nups" algorithm usually finds a good solution efficiently. See how it performs on this problem.

options = optimoptions("patternsearch",Algorithm="nups");
[x,fval,eflag,output] = patternsearch(fun,x0,A,b,[],[],[],[],[],options)
patternsearch stopped because the mesh size was less than options.MeshTolerance.

x = 1×2

    0.0088  -10.0088


fval = 
35.2488

eflag = 
1

output = struct with fields:
         function: @(t)sawtoothxy(t(1),t(2))
      problemtype: 'linearconstraints'
       pollmethod: 'nups'
    maxconstraint: 0
     searchmethod: []
       iterations: 61
        funccount: 180
         meshsize: 6.7893e-07
         rngstate: [1×1 struct]
          message: 'patternsearch stopped because the mesh size was less than options.MeshTolerance.'

The solver does not find the global minimum. However, choosing a different start point gives a better result.

x0 = [100,-50] + 10*randn(1,2); % New start point
[x,fval,eflag,output] = patternsearch(fun,x0,A,b,[],[],[],[],[],options)
patternsearch stopped because the mesh size was less than options.MeshTolerance.

x = 1×2
10-6 ×

    0.0247   -0.3333


fval = 
7.1681e-13

eflag = 
1

output = struct with fields:
         function: @(t)sawtoothxy(t(1),t(2))
      problemtype: 'linearconstraints'
       pollmethod: 'nups'
    maxconstraint: 0
     searchmethod: []
       iterations: 70
        funccount: 161
         meshsize: 9.6669e-07
         rngstate: [1×1 struct]
          message: 'patternsearch stopped because the mesh size was less than options.MeshTolerance.'

This time, the solver reaches the global minimum. But for a different problem, you might need to restart the solver many times.

Restart Solver from Random Points

To search for a global minimum, try starting patternsearch repeatedly from random points. If the problem has finite bounds on all variables, you can use the random points

x0 = lb + rand(size(lb)).*(ub - lb);

These points are uniformly generated at random within the bounds.

The current problem does not have explicit bounds. However, the linear constraints give a rectangular area that you can sample. More usefully, you can infer bounds on the variables: -15≤x≤200 and -105≤y≤110. Choose the inferred bounds as

lb = [-15,-105];
ub = [200,110];

Repeatedly choose random points using the uniform sampler to try to obtain the best result.

N = 20; % Number of attempts to make.
rng(600) % Set a point that does not lead to the global minimum.
x0 = lb + rand(size(lb)).*(ub - lb);
options.Display = "none"; % Suppress output.
[x,fval,eflag,output] = patternsearch(fun,x0,A,b,[],[],[],[],[],options);
disp(fval) % Starting function value.
   35.2488

for i = 2:N
    x0 = lb + rand(size(lb)).*(ub - lb);
    [x2,fval2,eflag2,output2] = patternsearch(fun,x0,A,b,[],[],[],[],[],options);
    if fval2 < fval % Copy results to x, fval, eflag, output
        x = x2;
        fval = fval2;
        eflag = eflag2;
        output = output2;
    end
end
disp(fval)
     0

The solver reaches the global solution. For other ways of searching for a global minimum, see Methods of Generating Start Points.

See Also

Topics