Simplex Search Source Code - description

# Simplex Search Source Code

A Simplex Search is a form of Direct Search, which is a class of non-derivative-based optimization methods. A simplex search uses a nondegenerate simplex as its design for function sampling.

The three Simplex Search classes were conceived as part of an undergraduate honors thesis by Adam Gurson under advisor Virginia Torczon in the Computer Science Department of the College of William and Mary. The project required access to a variety of simplex search algorithms that could be easily tested and compared. The classes implement three simplex search algorithms (discussed below): The simplex method of Spendley, Hext, and Himsworth; the algorithm by Nelder and Mead, as specified by Lagarias, Reeds, Wright, and Wright; and a sequential version of Torczon's multidimensional search. While each has its specific rules for simplex manipulation, they all follow roughly the same basic heuristic. Each begins with a nondegenerate simplex. The function is then evaluated at the vertices of the simplex. As the algorithms proceed, they keep track of the point at which the best function value is found. For each iteration of the algorithm, they attempt to replace simplex vertices that yield high function values with vertices whose new function values are lower. After each iteration, the vertices that replaced those in the starting simplex leave us with a new simplex. The goal of all three algorithms is to move the simplex, by replacing vertices, into the neighborhood of a minimizer.

The Vector and Matrix classes, built from the Template Numerical Toolkit, implement vectors and matrices in a fashion suitable for mathematical computation. The Template Numerical Toolkit was originally created by Roldan Pozo, who is a research scientist in the Mathematical and Computational Sciences Division of the Information Technology Laboratory at the National Institute of Standards and Technology (NIST). This version was modified by Chris Siefert from the College of William and Mary, to add increased functionality. New features include 2-norm computation (based on the BLAS implementation of drnm2) for vectors, alternative methods for accessing matrix data.

The SHHSearch class performs a simple search on trial points located in the positions described by the vertices of a rigid simplex. We adapted the classic Spendley, Hext, and Himsworth algorithm, originally designed for evolutionary operation rather than optimization, to better suit the purpose of numerical optimization. The underlying heuristic is to flip (by replacing vertices) the simplex along the gradient of the function, with the goal of locating a valid local minimizer. We use the notion of point age to decide when to shrink the simplex to narrow in on the exact location of a minimizer. For an in-depth discussion of the algorithm, see Adam Gurson's Honors Thesis.

The SMDSearch class is based on the algorithm first developed by Torczon for execution on a multi-processor machine. The version included here is adapted for single-processor machines. The main difference between the SMDSearch algorithm and the other simplex searches is that, during a reflection step, rather than simply reflecting the vertex with greatest function value through the centroid of the other vertices, the Multi-Directional Search instead reflects all vertices through the best vertex. As a result, each iteration of the SMDSearch actually involves two simplices: a primary simplex and its reflection through the current best minimizer. The algorithm uses a greedy approach to function sampling: function values are calculated only until a new minimizer is found. If a new minimizer is identified in the reflection simplex, the reflection simplex becomes the new primary simplex with the newly identified minimizer as its basepoint, and a new reflection simplex is constructed. If all 2n + 1 points in both primary and reflection simplices have been evaluated with no improvement found, the primary simplex shrinks around the basepoint, and a new iteration is begun, with a smaller simplex. More details are in Adam's Honors Thesis.

The objective files contain declarations and user-defined specifications including initial step length. objective.cc implements the functions declared in objective.h, including fcn() and Evalf(). The first of these, fcn(), is a wrapper for Evalf(), which uses the Unix fork-exec paradigm to call the actual function to be minimized. We have provided a sample function in the source file shekel.cc, so the user may become familiar with this process. Before running the test programs, the user must compile the shekel program on the same operating system in which the simplex searches will be compiled. Using g++, for example, the command would be:

g++ -Wall -g shekel.cc -o shekel

These test drivers offer examples of declaring the necessary classes, running the searches, and getting at some of the data possibilities. The searches run are simple, and the format is interactive. If users plan to use our shekel.cc file, they should first compile a binary executable as discussed above. Instructions for compiling each of the test programs can be found at the top of its source file.

An example of a run of NMtest is below (the '%' is just the prompt character):

% ./NMtest
SEED=23185
RUNS=3
dim=4
maxCalls=-1
Fname=shekel
Lower bound (space separated)=0 0 0 0
Upper bound (space separated)=9 9 9 9

NMSearch calls=225, Value=-10.5364, TolHit=1, Point=4.00076 4.00056 3.99963 3.99954
SHEKEL FAILED!!!!
NMSearch calls=152, Value=-2.87114, TolHit=1, Point=5.99906 5.99735 5.99828 5.99658
NMSearch calls=160, Value=-10.5364, TolHit=1, Point=4.00073 4.00053 3.99962 3.99956
INITIAL STEP LENGTH = 2

Note that for this function the user must use dim = 4, and should keep the lower and upper bounds within [0, 10]^4. The random seed must be a positive integer. The function evaluation must be printed to standard output. Any other output statements must be sent either to standard error or to a file. A value of -1 for maxCalls means that the search will terminate based only on the refinement of the design, rather than because it has reached its maximum allowable number of function calls.

Also note that the statement "SHEKEL FAILED" in the listing above does not indicate that the program is not working; it just means that the simplex has gone out of bounds at that point. The algorithm will flip (or otherwise adapt) the simplex to remedy this.