| Genetic Algorithm and Direct Search Toolbox™ | |
| On this page… |
|---|
Mesh Expansion and Contraction Setting Tolerances for the Solver |
Note All examples use the generalized pattern search (GPS) algorithm, but can be applied to the MADS algorithm as well. |
At each iteration, the pattern search polls the points in the current mesh—that is, it computes the objective function at the mesh points to see if there is one whose function value is less than the function value at the current point. How Pattern Search Works provides an example of polling. You can specify the pattern that defines the mesh by the Poll method option. The default pattern, GPS Positive basis 2N, consists of the following 2N directions, where N is the number of independent variables for the objective function.
[1 0 0...0]
[0 1 0...0]
...
[0 0 0...1]
[–1 0 0...0]
[0 –1 0...0]
[0 0 0...–1].
For example, if the objective function has three independent variables, the GPS Positive basis 2N, consists of the following six vectors.
[1 0 0]
[0 1 0]
[0 0 1]
[–1 0 0]
[0 –1 0]
[0
0 –1].
Alternatively, you can set Poll method to GPS Positive basis NP1, the pattern consisting of the following N + 1 directions.
[1 0 0...0]
[0 1 0...0]
...
[0 0 0...1]
[–1 –1
–1...–1].
For example, if objective function has three independent variables, the GPS Positive basis Np1, consists of the following four vectors.
[1 0 0]
[0 1 0]
[0 0 1]
[–1 –1 –1].
A pattern search will sometimes run faster using GPS Positive basis Np1 rather than the GPS Positive basis 2N as the Poll method, because the algorithm searches fewer points at each iteration. Although not being addressed in this example, the same is true when using the MADS Positive basis Np1 over the MADS Positive basis 2N. For example, if you run a pattern search on the example described in Example — A Linearly Constrained Problem, the algorithm performs 2080 function evaluations with GPS Positive basis 2N, the default Poll method, but only 1413 function evaluations using GPS Positive basis Np1.
However, if the objective function has many local minima, using GPS Positive basis 2N as the Poll method might avoid finding a local minimum that is not the global minimum, because the search explores more points around the current point at each iteration.
By default, if the pattern search finds a mesh point that improves the value of the objective function, it stops the poll and sets that point as the current point for the next iteration. When this occurs, some mesh points might not get polled. Some of these unpolled points might have an objective function value that is even lower than the first one the pattern search finds.
For problems in which there are several local minima, it is sometimes preferable to make the pattern search poll all the mesh points at each iteration and choose the one with the best objective function value. This enables the pattern search to explore more points at each iteration and thereby potentially avoid a local minimum that is not the global minimum. You can make the pattern search poll the entire mesh setting Complete poll to On in Poll options.
As an example, consider the following function.
The following figure shows a plot of the function.
The global minimum of the function occurs at (0, 0), where its value is -25. However, the function also has a local minimum at (0, 9), where its value is -16.
To create an M-file that computes the function, copy and paste the following code into a new M-file in the MATLAB Editor.
function z = poll_example(x)
if x(1)^2 + x(2)^2 <= 25
z = x(1)^2 + x(2)^2 - 25;
elseif x(1)^2 + (x(2) - 9)^2 <= 16
z = x(1)^2 + (x(2) - 9)^2 - 16;
else z = 0;
end
Then save the file as poll_example.m in a directory on the MATLAB path.
To run a pattern search on the function, enter the following in the Optimization Tool:
Set Solver to patternsearch.
Set Objective function to @poll_example.
Set Start point to [0 5].
Set Level of display to Iterative in the Display to command window options.
Click Start to run the pattern search with Complete poll set to Off, its default value. The Optimization Tool displays the results in the Run solver and view results pane, as shown in the following figure.
The pattern search returns the local minimum at (0, 9). At the initial point, (0, 5), the objective function value is 0. At the first iteration, the search polls the following mesh points.
f((0, 5) + (1, 0)) = f(1, 5) = 0
f((0, 5) + (0, 1)) = f(0, 6) = -7
As soon as the search polls the mesh point (0, 6), at which the objective function value is less than at the initial point, it stops polling the current mesh and sets the current point at the next iteration to (0, 6). Consequently, the search moves toward the local minimum at (0, 9) at the first iteration. You see this by looking at the first two lines of the command line display.
Iter f-count MeshSize f(x) Method
0 1 1 0 Start iterations
1 3 2 -7 Successful Poll
Note that the pattern search performs only two evaluations of the objective function at the first iteration, increasing the total function count from 1 to 3.
Next, set Complete poll to On and click Start. The Run solver and view results pane displays the following results.
This time, the pattern search finds the global minimum at (0, 0). The difference between this run and the previous one is that with Complete poll set to On, at the first iteration the pattern search polls all four mesh points.
f((0, 5) + (1, 0)) = f(1, 5) = 0
f((0, 5) + (0, 1)) = f(0, 6) = -6
f((0, 5) + (-1, 0)) = f(-1, 5) = 0
f((0, 5) + (0, -1)) = f(0, 4) = -9
Because the last mesh point has the lowest objective function value, the pattern search selects it as the current point at the next iteration. The first two lines of the command-line display show this.
Iter f-count MeshSize f(x) Method
0 1 1 0 Start iterations
1 5 2 -9 Successful Poll
In this case, the objective function is evaluated four times at the first iteration. As a result, the pattern search moves toward the global minimum at (0, 0).
The following figure compares the sequence of points returned when Complete poll is set to Off with the sequence when Complete poll is On.
In addition to polling the mesh points, the pattern search algorithm can perform an optional step at every iteration, called search. At each iteration, the search step applies another optimization method to the current point. If this search does not improve the current point, the poll step is performed.
The following example illustrates the use of a search method on the problem described in Example — A Linearly Constrained Problem. To set up the example, enter the following commands at the MATLAB prompt to define the initial point and constraints.
x0 = [2 1 0 9 1 0]; Aineq = [-8 7 3 -4 9 0 ]; bineq = [7]; Aeq = [7 1 8 3 3 3; 5 0 5 1 5 8; 2 6 7 1 1 8; 1 0 0 0 0 0]; beq = [84 62 65 1];
Then enter the settings shown in the following figures in the Optimization Tool.
For comparison, click Start to run the example without a search method. This displays the plots shown in the following figure.
To see the effect of using a search method, select GPS Positive Basis Np1 in the Search method field in Search options. This sets the search method to be a pattern search using the pattern for GPS Positive basis Np1. Then click Start to run the genetic algorithm. This displays the following plots.
Note that using the search method reduces the total function evaluations by almost 50 percent—from 1889 to 981—and reduces the number of iterations from 250 to 78.
The Expansion factor and Contraction factor options, in Mesh options, control how much the mesh size is expanded or contracted at each iteration. With the default Expansion factor value of 2, the pattern search multiplies the mesh size by 2 after each successful poll. With the default Contraction factor value of 0.5, the pattern search multiplies the mesh size by 0.5 after each unsuccessful poll.
You can view the expansion and contraction of the mesh size during the pattern search by selecting Mesh size in the Plot functions pane.
To also display the values of the mesh size and objective function at the command line, set Level of display to Iterative in the Display to command window options.
When you run the example described in Example — A Linearly Constrained Problem, the Optimization Tool displays the following plot.
To see the changes in mesh size more clearly, change the y-axis to logarithmic scaling as follows:
Updating these settings in the MATLAB Property Editor will show the plot in the following figure.
The first 37 iterations result in successful polls, so the mesh sizes increase steadily during this time. You can see that the first unsuccessful poll occurs at iteration 38 by looking at the command-line display for that iteration.
36 39 6.872e+010 3486 Successful Poll 37 40 1.374e+011 3486 Successful Poll 38 43 6.872e+010 3486 Refine Mesh
Note that at iteration 37, which is successful, the mesh size doubles for the next iteration. But at iteration 38, which is unsuccessful, the mesh size is multiplied 0.5.
To see how Expansion factor and Contraction factor affect the pattern search, make the following changes:
Set Expansion factor to 3.0.
Set Contraction factor to 0.75.
Then click Start. The Run solver and view results pane shows that the final point is approximately the same as with the default settings of Expansion factor and Contraction factor, but that the pattern search takes longer to reach that point.
The algorithm halts because it exceeds the maximum number of iterations, whose value you can set in the Max iteration field in the Stopping criteria options. The default value is 100 times the number of variables for the objective function, which is 6 in this example.
When you change the scaling of the y-axis to logarithmic, the mesh size plot appears as shown in the following figure.
Note that the mesh size increases faster with Expansion factor set to 3.0, as compared with the default value of 2.0, and decreases more slowly with Contraction factor set to 0.75, as compared with the default value of 0.5.
The mesh accelerator can make a pattern search converge faster to an optimal point by reducing the number of iterations required to reach the mesh tolerance. When the mesh size is below a certain value, the pattern search contracts the mesh size by a factor smaller than the Contraction factor factor.
Note It is recommended to only use the mesh accelerator for problems in which the objective function is not too steep near the optimal point, or you might lose some accuracy. For differentiable problems, this means that the absolute value of the derivative is not too large near the solution. |
To use the mesh accelerator, set Accelerator to On in the Mesh options. When you run the example described in Example — A Linearly Constrained Problem, the number of iterations required to reach the mesh tolerance is 246, as compared with 270 when Accelerator is set to Off.
You can see the effect of the mesh accelerator by setting Level of display to Iterative in Display to command window. Run the example with Accelerator set to On, and then run it again with Accelerator set to Off. The mesh sizes are the same until iteration 226, but differ at iteration 227. The MATLAB Command Window displays the following lines for iterations 226 and 227 with Accelerator set to Off.
Iter f-count MeshSize f(x) Method 226 1501 6.104e-005 2189 Refine Mesh 227 1516 3.052e-005 2189 Refine Mesh
Note that the mesh size is multiplied by 0.5, the default value of Contraction factor.
For comparison, the Command Window displays the following lines for the same iteration numbers with Accelerator set to On.
Iter f-count MeshSize f(x) Method 226 1501 6.104e-005 2189 Refine Mesh 227 1516 1.526e-005 2189 Refine Mesh
In this case the mesh size is multiplied by 0.25.
Typically, at any given iteration of a pattern search, some of the mesh points might coincide with mesh points at previous iterations. By default, the pattern search recomputes the objective function at these mesh points even though it has already computed their values and found that they are not optimal. If computing the objective function takes a long time—say, several minutes—this can make the pattern search run significantly longer.
You can eliminate these redundant computations by using a cache, that is, by storing a history of the points that the pattern search has already visited. To do so, set Cache to On in Cache options. At each poll, the pattern search checks to see whether the current mesh point is within a specified tolerance, Tolerance, of a point in the cache. If so, the search does not compute the objective function for that point, but uses the cached function value and moves on to the next point.
Note When Cache is set to On, the pattern search might fail to identify a point in the current mesh that improves the objective function because it is within the specified tolerance of a point in the cache. As a result, the pattern search might run for more iterations with Cache set to On than with Cache set to Off. It is generally a good idea to keep the value of Tolerance very small, especially for highly nonlinear objective functions. |
To illustrate this, select Best function value and Function count in the Plot functions pane and run the example described in Example — A Linearly Constrained Problem with Cache set to Off. After the pattern search finishes, the plots appear as shown in the following figure.
Note that the total function count is 2080.
Now, set Cache to On and run the example again. This time, the plots appear as shown in the following figure.
This time, the total function count is reduced to 1973.
Tolerance refers to how small a parameter, such a mesh size, can become before the search is halted or changed in some way. You can specify the value of the following tolerances:
Mesh tolerance — When the current mesh size is less than the value of Mesh tolerance, the algorithm halts.
X tolerance — After a successful poll, if the distance from the previous best point to the current best point is less than the value of X tolerance, the algorithm halts.
Function tolerance — After a successful poll, if the difference between the function value at the previous best point and function value at the current best point is less than the value of Function tolerance, the algorithm halts.
Nonlinear constraint tolerance — The algorithm treats a point to be feasible if constraint violation is less than TolCon.
Bind tolerance — Bind tolerance applies to constrained problems and specifies how close a point must get to the boundary of the feasible region before a linear constraint is considered to be active. When a linear constraint is active, the pattern search polls points in directions parallel to the linear constraint boundary as well as the mesh points.
Usually, you should set Bind tolerance to be at least as large as the maximum of Mesh tolerance, X tolerance, and Function tolerance.
The following example illustrates of how Bind tolerance affects a pattern search. The example finds the minimum of
subject to the constraints
–11x1 +
10x2 ≤ 10
10x1 –
10x2 ≤ 10.
Note that you can compute the objective function using the function norm. The feasible region for the problem lies between the two lines in the following figure.
To run the example, enter optimtool and choose patternsearch to open the Optimization Tool. Enter the following information:
Set Objective function to @(x) norm(x).
Set Start point to [-1.001 -1.1].
Select Mesh size in the Plot functions pane.
Set Level of display to Iterative in the Display to command window options.
Then click Start to run the pattern search.
The display in the MATLAB Command Window shows that the first four polls are unsuccessful, because the mesh points do not lie in the feasible region.
Iter f-count MeshSize f(x) Method
0 1 1 1.487 Start iterations
1 1 0.5 1.487 Refine Mesh
2 1 0.25 1.487 Refine Mesh
3 1 0.125 1.487 Refine Mesh
4 1 0.0625 1.487 Refine Mesh
The pattern search contracts the mesh at each iteration until one of the mesh points lies in the feasible region. The following figure shows a close-up of the initial point and mesh points at iteration 5.
The top mesh point, which is (-1.001, -1.0375), has a smaller objective function value than the initial point, so the poll is successful.
Because the distance from the initial point to lower boundary line is less than the default value of Bind tolerance, which is 0.0001, the pattern search does not consider the linear constraint 10x1 – 10x2 ≤ 10 to be active, so it does not search points in a direction parallel to the boundary line.
To see the effect of bind tolerance, change Bind tolerance to 0.01 and run the pattern search again.
This time, the display in the MATLAB Command Window shows that the first two iterations are successful.
Iter f-count MeshSize f(x) Method
0 1 1 1.487 Start iterations
1 2 2 0.7817 Successful Poll
2 3 4 0.6395 Successful Poll
Because the distance from the initial point to the boundary is less than Bind tolerance, the second linear constraint is active. In this case, the pattern search polls points in directions parallel to the boundary line 10x1 – 10x2 ≤ 10, resulting in successful poll. The following figure shows the initial point with two addition search points in directions parallel to the boundary.
The following figure compares the sequences of points during the first 20 iterations of the pattern search for both settings of Bind tolerance.
Note that when Bind tolerance is set to .01, the points move toward the optimal point more quickly. The pattern search requires only 90 iterations. When Bind tolerance is set to .0001, the search requires 124 iterations. However, when the feasible region does not contain very acute angles, as it does in this example, increasing Bind tolerance can increase the number of iterations required, because the pattern search tends to poll more points.
Suppose you want to minimize the simple objective function of two variables x1 and x2,
subject to the following nonlinear inequality constraints and bounds
Begin by creating the objective and constraint functions. First, create an M-file named simple_objective.m as follows:
function y = simple_objective(x) y = (4 - 2.1*x(1)^2 + x(1)^4/3)*x(1)^2 + x(1)*x(2) + (-4 + 4*x(2)^2)*x(2)^2;
The pattern search solver assumes the objective function will take one input x where x has as many elements as number of variables in the problem. The objective function computes the value of the function and returns that scalar value in its one return argument y.
Then create an M-file named simple_constraint.m containing the constraints:
function [c, ceq] = simple_constraint(x) c = [1.5 + x(1)*x(2) + x(1) - x(2); -x(1)*x(2) + 10]; ceq = [];
The pattern search solver assumes the constraint function will take one input x, where x has as many elements as the number of variables in the problem. The constraint function computes the values of all the inequality and equality constraints and returns two vectors, c and ceq, respectively.
Next, to minimize the objective function using the patternsearch function, you need to pass in a function handle to the objective function as well as specifying a start point as the second argument. Lower and upper bounds are provided as LB and UB respectively. In addition, you also need to pass a function handle to the nonlinear constraint function.
ObjectiveFunction = @simple_objective;
X0 = [0 0]; % Starting point
LB = [0 0]; % Lower bound
UB = [1 13]; % Upper bound
ConstraintFunction = @simple_constraint;
[x,fval] = patternsearch(ObjectiveFunction,X0,[],[],[],[],...
LB,UB,ConstraintFunction)
Optimization terminated: mesh size less than options.TolMesh
and constraint violation is less than options.TolCon.
x =
0.8122 12.3122
fval =
9.1324e+004Next, plot the results. Create an options structure using psoptimset that selects two plot functions. The first plot function psplotbestf plots the best objective function value at every iteration. The second plot function psplotmaxconstr plots the maximum constraint violation at every iteration.
Note You can also visualize the progress of the algorithm by displaying information to the Command Window using the 'Display' option. |
options = psoptimset('PlotFcns',{@psplotbestf,@psplotmaxconstr},'Display','iter');
[x,fval] = patternsearch(ObjectiveFunction,X0,[],[],[],[],LB,UB,ConstraintFunction,options)
max
Iter f-count f(x) constraint MeshSize Method
0 1 0 10 0.8919
1 5 113580 0 0.001 Increase penalty
2 24 91324.4 0 1e-005 Increase penalty
3 48 91324 0 1e-007 Increase penalty
Optimization terminated: mesh size less than options.TolMesh
and constraint violation is less than options.TolCon.
x =
0.8122 12.3122
fval =
9.1324e+004Best Objective Function Value and Maximum Constraint Violation at Each Iteration
Direct search often runs faster if you vectorize the objective and nonlinear constraint functions. This means your functions evaluate all the points in a poll or search pattern at once, with one function call, without having to loop through the points one at a time. Therefore, the option Vectorize = 'on' works only when CompletePoll or CompleteSearch are also set to 'on'.
If there are nonlinear constraints, the objective function and the nonlinear constraints all need to be vectorized in order for the algorithm to compute in a vectorized manner.
A vectorized objective function accepts a matrix as input and generates a vector of function values, where each function value corresponds to one row or column of the input matrix. patternsearch resolves the ambiguity in whether the rows or columns of the matrix represent the points of a pattern as follows. Suppose the input matrix has m rows and n columns:
If the initial point x0 is a column vector of size m, the objective function takes each column of the matrix as a point in the pattern and returns a vector of size n.
If the initial point x0 is a row vector of size n, the objective function takes each row of the matrix as a point in the pattern and returns a vector of size m.
If the initial point x0 is a scalar, the matrix has one row (m = 1, the matrix is a vector), and each entry of the matrix represents one point to evaluate.
Pictorially, the matrix and calculation are represented by the following figure.
Structure of Vectorized Functions
For example, suppose the objective function is
If the initial vector x0 is a column vector, such as [0;0], an M-file for vectorized evaluation is
function f = vectorizedc(x)
f = x(1,:).^4+x(2,:).^4-4*x(1,:).^2-2*x(2,:).^2...
+3*x(1,:)-.5*x(2,:);If the initial vector x0 is a row vector, such as [0,0], an M-file for vectorized evaluation is
function f = vectorizedr(x)
f = x(:,1).^4+x(:,2).^4-4*x(:,1).^2-2*x(:,2).^2...
+3*x(:,1)-.5*x(:,2);If you want to use the same objective (fitness) function for both pattern search and genetic algorithm, write your function to have the points represented by row vectors, and write x0 as a row vector. The genetic algorithm always takes individuals as the rows of a matrix. This was a design decision—the genetic algorithm does not require a user-supplied population, so needs to have a default format.
To minimize vectorizedc, enter the following commands:
options=psoptimset('Vectorized','on','CompletePoll','on');
x0=[0;0];
[x fval]=patternsearch(@vectorizedc,x0,...
[],[],[],[],[],[],[],options)MATLAB returns the following output:
Optimization terminated: mesh size less than options.TolMesh.
x =
-1.5737
1.0575
fval =
-10.0088Only nonlinear constraints need to be vectorized; bounds and linear constraints are handled automatically. If there are nonlinear constraints, the objective function and the nonlinear constraints all need to be vectorized in order for the algorithm to compute in a vectorized manner.
The same considerations hold for constraint functions as for objective functions: the initial point x0 determines the type of points (row or column vectors) in the poll or search. If the initial point is a row vector of size k, the matrix x passed to the constraint function has k columns. Similarly, if the initial point is a column vector of size k, the matrix of poll or search points has k rows. The figure Structure of Vectorized Functions may make this clear.
Your nonlinear constraint function returns two matrices, one for inequality constraints, and one for equality constraints. Suppose there are nc nonlinear inequality constraints and nceq nonlinear equality constraints. For row vector x0, the constraint matrices have nc and nceq columns respectively, and the number of rows is the same as in the input matrix. Similarly, for a column vector x0, the constraint matrices have nc and nceq rows respectively, and the number of columns is the same as in the input matrix. In figure Structure of Vectorized Functions, "Results" includes both nc and nceq.
Suppose that the nonlinear constraints are
Write an M-file for these constraints for row-form x0 as follows:
function [c ceq] = ellipsecosh(x) c(:,1)=x(:,1).^2/9+x(:,2).^2/4-1; c(:,2)=cosh(x(:,1))-x(:,2)-1; ceq=[];
Minimize vectorizedr (defined in Vectorized Objective Function) subject to the constraints ellipsecosh:
x0=[0,0];
options=psoptimset('Vectorized','on','CompletePoll','on');
[x fval]=patternsearch(@vectorizedr,x0,...
[],[],[],[],[],[],@ellipsecosh,options);MATLAB returns the following output:
Optimization terminated: mesh size less than options.TolMesh and constraint violation is less than options.TolCon. x = -1.3516 1.0612 fval = -9.5394
| Performing a Pattern Search from the Command Line | Parallel Computing with Pattern Search | |
| © 1984-2008- The MathWorks, Inc. - Site Help - Patents - Trademarks - Privacy Policy - Preventing Piracy - RSS |