Documentation |
Finding a start point in the basin of attraction of the global minimum can be difficult when the basin is small or when you are unsure of the location of the minimum. To solve this type of problem you can:
Add sensible bounds
Take a huge number of random start points
Make a methodical grid of start points
For an unconstrained problem, take widely dispersed random start points
This example shows these methods and some variants.
The function –sech(x) is nearly 0 for all |x| > 5, and –sech(0) = –1. The example is a two-dimensional version of the sech function, with one minimum at [1,1], the other at [1e5,-1e5]:
f(x,y) = –10sech(|x – (1,1)|) – 20sech(.0003(|x – (1e5,–1e5)|) – 1.
f has a global minimum of –21 at (1e5,–1e5), and a local minimum of –11 at (1,1).
Code for Generating the Figure
The minimum at (1e5,–1e5) shows as a narrow spike. The minimum at (1,1) does not show since it is too narrow.
The following sections show various methods of searching for the global minimum. Some of the methods are not successful on this problem. Nevertheless, you might find each method useful for different problems.
GlobalSearch and MultiStart cannot find the global minimum using default global options, since the default start point components are in the range (–9999,10001) for GlobalSearch and (–1000,1000) for MultiStart.
With additional bounds of –1e6 and 1e6 in problem, GlobalSearch usually does not find the global minimum:
x1 = [1;1];x2 = [1e5;-1e5]; f = @(x)-10*sech(norm(x(:)-x1)) -20*sech((norm(x(:)-x2))*3e-4) -1; opts = optimoptions(@fmincon,'Algorithm','active-set'); problem = createOptimProblem('fmincon','x0',[0,0],'objective',f,... 'lb',[-1e6;-1e6],'ub',[1e6;1e6],'options',opts); gs = GlobalSearch; rng(14,'twister') % for reproducibility [xfinal fval] = run(gs,problem) GlobalSearch stopped because it analyzed all the trial points. All 32 local solver runs converged with a positive local solver exit flag. xfinal = 1.0000 1.0000 fval = -11.0000
To find the global minimum, you can search more points. This example uses 1e5 start points, and a MaxTime of 300 s:
gs.NumTrialPoints = 1e5; gs.MaxTime = 300; [xg fvalg] = run(gs,problem) GlobalSearch stopped because maximum time is exceeded. GlobalSearch called the local solver 2186 times before exceeding the clock time limit (MaxTime = 300 seconds). 1943 local solver runs converged with a positive local solver exit flag. xg = 1.0e+04 * 10.0000 -10.0000 fvalg = -21.0000
In this case, GlobalSearch found the global minimum.
Alternatively, you can search using MultiStart with many start points. This example uses 1e5 start points, and a MaxTime of 300 s:
ms = MultiStart(gs); [xm fvalm] = run(ms,problem,1e5) MultiStart stopped because maximum time was exceeded. MultiStart called the local solver 17266 times before exceeding the clock time limit (MaxTime = 300 seconds). 17266 local solver runs converged with a positive local solver exit flag. xm = 1.0000 1.0000 fvalm = -11.0000
In this case, MultiStart failed to find the global minimum.
You can also use MultiStart to search an unbounded region to find the global minimum. Again, you need many start points to have a good chance of finding the global minimum.
The first five lines of code generate 10,000 widely dispersed random start points using the method described in Widely Dispersed Points for Unconstrained Components. newprob is a problem structure using the fminunc local solver and no bounds:
rng(0,'twister') % for reproducibility u = rand(1e4,1); u = 1./u; u = exp(u) - exp(1); s = rand(1e4,1)*2*pi; stpts = [u.*cos(s),u.*sin(s)]; startpts = CustomStartPointSet(stpts); opts = optimoptions(@fminunc,'Algorithm','quasi-newton'); newprob = createOptimProblem('fminunc','x0',[0;0],'objective',f,... 'options',opts); [xcust fcust] = run(ms,newprob,startpts) MultiStart completed the runs from all start points. All 10000 local solver runs converged with a positive local solver exit flag. xcust = 1.0e+05 * 1.0000 -1.0000 fcust = -21.0000
In this case, MultiStart found the global minimum.
You can also use a grid of start points instead of random start points. To learn how to construct a regular grid for more dimensions, or one that has small perturbations, see Uniform Grid or Perturbed Grid.
xx = -1e6:1e4:1e6; [xxx yyy] = meshgrid(xx,xx); z = [xxx(:),yyy(:)]; bigstart = CustomStartPointSet(z); [xgrid fgrid] = run(ms,newprob,bigstart) MultiStart completed the runs from all start points. All 10000 local solver runs converged with a positive local solver exit flag. xgrid = 1.0e+004 * 10.0000 -10.0000 fgrid = -21.0000
In this case, MultiStart found the global minimum.
Making a regular grid of start points, especially in high dimensions, can use an inordinate amount of memory or time. You can filter the start points to run only those with small objective function value.
To perform this filtering most efficiently, write your objective function in a vectorized fashion. For information, see Write a Vectorized Function or Vectorize the Objective and Constraint Functions. The following function handle computes a vector of objectives based on an input matrix whose rows represent start points:
x1 = [1;1];x2 = [1e5;-1e5]; g = @(x) -10*sech(sqrt((x(:,1)-x1(1)).^2 + (x(:,2)-x1(2)).^2)) ... -20*sech(sqrt((x(:,1)-x2(1)).^2 + (x(:,2)-x2(2)).^2))-1;
Suppose you want to run the local solver only for points where the value is less than –2. Start with a denser grid than in MultiStart with a Regular Grid of Start Points, then filter out all the points with high function value:
xx = -1e6:1e3:1e6; [xxx yyy] = meshgrid(xx,xx); z = [xxx(:),yyy(:)]; idx = g(z) < -2; % index of promising start points zz = z(idx,:); smallstartset = CustomStartPointSet(zz); opts = optimoptions(@fminunc,'Algorithm','quasi-newton','Display','off'); newprobg = createOptimProblem('fminunc','x0',[0,0],... 'objective',g,'options',opts); % row vector x0 since g expects rows [xfew ffew] = run(ms,newprobg,smallstartset) MultiStart completed the runs from all start points. All 2 local solver runs converged with a positive local solver exit flag. xfew = 100000 -100000 ffew = -21
In this case, MultiStart found the global minimum. There are only two start points in smallstartset, one of which is the global minimum.