User's Manual for DirectSearch and its Derived Classes

P.L. (Anne) Shepherd


Contents

System requirements

This software was designed using Linux (Red Hat 6.2) and tested using various flavors of UNIX--specifically, Linux, SunOS, and AIX. Our usual compiler was some version of gcc, mostly gcc version egcs-2.91.66 19990314/Linux (egcs-1.1.2 release). If you plan to use the EvalF() function to call your objective function, you will have to use some form of UNIX. If you are using a non-UNIX operating system, it is likely that you will have to make modifications to this software.

Where to find the software

You can find the DirectSearch software at QQQQQQQ. Just follow the directions given there.

Installation

Once you have the gzipped tar file in the directory where you want it to be, Unzip and untar the file. You can do this with the command: 

tar zxvf DirectSearch.tar.gz

If your version of tar doesn't support the 'z' option (gtar does), then do the following: 

gunzip DirectSearch.tar.gz (This will create DirectSearch.tar) 

tar xvf DirectSearch.tar

A directory called direct_code, will be created in the directory in which you unzip the tar file. Type 

cd ./direct_code
to get to the source code.

A Makefile is included with the software. Typing ``make'' at the command prompt in the directory where the source code is located should make the demos main_eval, main_default, main_arb, main_fail, main_hy, main_obj, and main_spec. The source code files for these main programs include compilation instructions in the header comments in case you wish to compile any of them by hand.

Using allocator templates found in Dyn_alloc.h

If your compiler doesn not fully support ANSI Standard C++, you may get error messages when you try to compile this software.

If this happens, recompile with the flag -DDOLD_ALLOC. This will tell the compiler to use the old-style allocators, without try-blocks. If you are using the Makefile, comment out the CFLAGS line that does not have this flag, then comment in (by removing the pound sign) the line below it, which includes the -DDOLD_ALLOC flag.

Hooking your objective function up to our searches

There are three ways of doing this: wrap your function directly in a file called objective.cc; send the function name in as a parameter to a constructor (but some rewriting will almost certainly have to be done there as well); or use our EvalF() function to call the objective directly. The files objective.h and objective.cc contain functions that can be used for the various options.


Wrapping your function


The program main_default.cc shows how to use one of our searches using an objective that has been wrapped by a function called fcn(). Since ``fcn'' is the function name used as a default in all our searches, the short constructor can be used in this case. fcn() is defined in objective.cc as follows: 

void fcn(long vars, Vector<double> &x, double & f, 
         bool & flag, void* nothing)
{
     
  if (vars < 1) {
      cerr << "\nError: Dimension cannot be less than 1!!\n";
      exit(1);
  }
  f = mySin(x);
  flag = true;
  return;
}

Note that all it does is to call the function mySin(), defined later in objective.cc. It is easy enough to change fcn() to wrap any function that takes in a Vector and returns a function value. A function that uses an array instead of a Vector can be easily converted, using methods included in vec.h For example, thisVector.begin() points to the beginning of the array that represents the vector. An array can be made into a Vector as in the following example: 

. . .

int n = 4;         // dimension of the search

// It doesn't matter how you get your array, but 
// we'll do it like this:
double startVals[] =  {3, -1, 0, 1}; 


// We'll initialize an n-entry Vector whose value is startVal,
// and use it as our starting point.
Vector<double>minVec(n, startVals);

// Now we can call the constructor.
HJSearch HJS(n, minVec);

. . .


Making the function name an argument to the constructor.

The program main_spec.cc shows how to send the function name in as an argument to the constructor. In this example, we optimize Powell's quartic function (called powell() in our objective.cc file) , with initial settings as used by Nelder and Mead.

The call to the constructor is: 

 
NMSearch NM(n, minVec, 0.5, 1.0, 0.5, 2.0, startstep,
            endstep, powell, NULL);

where the third through sixth parameters are actually the same as the default settings for the associated fields.

In objective.cc, the function powell() is formulated so that it takes arguments of the same type and in the same order as in the DirectSearch field void(*fcn_name)(long dim, Vector$<$double$>$ &x, double & function, bool & success, void* an_object).

It is also possible to optimize a function from a separate C++ class, using the same constructor, with the last parameter non-NULL, as follows:

In the case of a CompassSearch, for example, the call to the constructor might look more like this: 

 
  ...

  // construct an instance of the class that holds the objective
  // function.
  MyClass* MyClassObject  = new MyClass();
  
  CompassSearch CMPSearcher(n, minVec, start, tol, 
                            MyFriendFunction, 
                            (void *)MyClassObject);

In this example, the field fcn_name points to MyFriendFunction, a ``friend'' function of the class MyClass. MyFriendFunction wraps the actual objective function, a member function of MyClass. The wrapping is necessary to get around the fact that C++ does not allow passing of pointers to class methods without the class object. So we include as parameters both a pointer to the class object and a pointer to the friend function that manipulates it.

I have provided a toy class trigclass and an example program main_obj.cc to illustrate this process.


Using fcn_eval() and EvalF() to call a binary executable


If the objective you are optimizing can be adapted so it takes a point's entries on its standard input and reports the function value on its standard output, you can use this option. Here, you can write your function in a language other than C/C++ if you wish, as long as the executable has the behavior just described.

The function fcn_eval() (in our objective.cc file) calls EvalF(), uses the UNIX fork and exec calls to spawn a child process and superimpose the user's objective function upon the child process. The parent and child processes communicate via pipes. The parent process (whose flow of control is within EvalF() ) sends the child process (the user's executable) the point it wants evaluated by printing it to its (the parent's) standard output, which has been redirected to a pipe. Your program (which will be the child process) must do the evaluation and then output the function value on its standard output.

Be advised that if your program prints out anything else to standard output, the whole thing will probably crash. Therefore, if you need or want to print anything out from your objective function, be sure to print either to standard error or to a file.

If the function evaluation is successful, EvalF() expects a return value of zero; if unsuccessful, a nonzero value should be returned.

Our example program main_eval.cc shows how to set up a search using this option. If you are using the Makefile and wigh to compile only what is needed for this program, type ``make eval'' at the command prompt. This will make the main_eval and standalone executables. If compiling by hand, you must first compile the program standalone.cc and name the executable standalone. To do this using g++, enter the following at the command line: 

g++ -g -Wall standalone.cc -o standalone

Customizing the Searches

Many of the fields of the searches can be mnaipulated by the user. See Appendix A for a list of all the fields and the accessors and mutators associated with them.

Scaling the searches

Throughout these searches we use a default of 2.0 as the starting length of an edge of the design. Depending on the scale of the problem being optimized, the user will wish to rescale the search. To do this, use the method void SetInitialStepLength(double & stepLen) for the Pattern searches and void SetStartingEdgeLengths(const Vector<double> & lengths) for the Simplex searches.

Telling the searches when to stop

Using a function call budget

By default, all the searches will terminate when their step size (or simplex edge length) falls below 10e-8, approximately the square root of machine epsilon. But there are times when the user may want to specify a maximum number of allowable calls (a ``function call budget.'') for a particular search. This is useful, for example, in comparing the performance of different searches in different situations; and may be necessary where function calls are very expensive.

To manage the call budget, all searches contain the fields maxCalls and exact_count, along with their associated accessors and mutators.

The field maxCalls is the maximum number of calls allowed. This is set by default to NO_MAX (which is -1), but can be set by the user. GetMaxCalls() returns the value of maxCalls. The user can define the call budget with the mutators SetMaxCalls(longcalls) and SetMaxCallsExact(long calls).

The boolean field exact_count indicates that the call budget is firm. This is set by default to ``true'' in DirectSearch constructors, and reset to ``false'' in all the Pattern Searches (including SMDSearch, which is a member of SimplexSearch class). If exact_count is true and a budget has been set (i.e. maxCalls does not equal NO_MAX), the search will terminate immediately as soon as functionCalls is equal to maxCalls. If there is a set budget and exact_count is false, however, the search will terminate only after a redunction in delta (a shrinking step) has taken place. The rationale is that, in a Pattern Search, a reduction in delta occurs only when a search of the points in the pattern yields no improvement, so this is a natural measure of progress to test for convergence. [1] There is no use in halting the search before such a reduction unless very strict call budgeting is necessary.

In the case of the Nelder-Mead search, we strongly recommend against setting exact_count to false. The nature of the Nelder-Mead search is such that it does not necessarily ever make a shrinking step. The simplex in this search deforms to suit the search. [2,5] It can--and often does--converge without ever shrinking. (But see the section below on stopping criteria in the Simplex Searches.)

The methods SetExact() and SetInexact() set the value of exact_count to true and false, respectively; and the boolean method IsExact() returns the value of exact_count. SetMaxCalls(long calls) sets maxCalls to the value of calls and at the same time sets exact_count to true.

Stopping criteria in the Simplex Searches

Although pattern searches are usually halted based on the refinement of the underlying lattice (presuming the search is within its call budget), this is not necessarily the case for simplex searches. The stopping criterion for the search originally described by Nelder and Mead [5] is a standard-deviation test based on the values of the simplex points. For more information, see [2,5].

There are times, however, when it is preferable to stop the search based rather upon refinement of the simplex; in other words, an equivalent to the ``delta'' test used in the pattern searches.

To manage this choice, we have included the boolean field Stop_on_std to indicate which criterion is being used by a particular simplex search. Another field, delta, measures the refinement of the simplex. The method void Set_Stop_on_std() sets Stop_on_std to true; void Set_Stop_on_delta() sets it to false. The boolean accessor Is_Stop_on_std() const returns its value;


In the case of SMDSearch class, delta means precisely what it does in the case of any other pattern search. For the SHHSearch class, which used a simplex that does not deform, delta has the same meaning as in the pattern searches and is calculated in the same way: it is set as the length of one side of the simplex at initialization; when the simplex shrinks, delta is reduced by the same scale. The default setting for Stop_on_std for SHHSearch and SHHSearch is ''false.''

For NMSearch class, however, the default is ``true.'' We made this decision based, first, on the description of the algorithm by Nelder and Mead [5], and, second, because a ``delta'' test cannot be as easily done in the NMSearch. Because the simplex in this search deforms to adapt to the search conditions, it becomes expensive to calculate NMdelta, the NMSearch equivalent of ``delta.'' The cost of each calculation is O($N^3$), where nis the dimension of the search. The standard-deviation test, by comparison, is O($N$).

If Stop_on_std is set to false in an NMSearch, the standard-deviation test is used until that test is met. Then (and in subsequent iterations of the same search) the more expensive ``delta'' test will be used. We realize that there will be times when the standard deviation test will actually take more function calls than the ``delta'' test, but this is not typically the case, and we make this compromise for the sake of efficiency.

For purposes of uniformity, SMDSearch and SHHSearch classes have the option of using either stopping criterion as well. NMSearch, however, is the only search that has Stop_on_std set to true as a default.


Please be advised that using the standard-deviation test can on rare occasions result in premature convergence. An example of this (using an SHHSearch) can be seen by compiling and running the program main_fail.cc. In this example, the search becomes stuck on a level set of the objective function. That is why we recommend using the ``delta'' stopping test whenever possible.


Getting around the problem by using Hybrid_NMSearch class


Another approach is to use a Hybrid_NMSearch. This class was written specifically to get around the high cost of calculating NMdelta. In this search, an NMSearch is constructed and run until the standard-deviation test is met (Stop_on_std is set to true in all constructors). Then, if there are function calls left in the budget and the ``delta'' test has not yet been satisfied, an EdHJSearch object is constructed using settings taken from the end state of the NMSearch.

This obviates the need for the more expensive calculation required by the NMSearch. Be advised, though, that any Nelder-Mead search has a chance of converging early (and falsely), due to the possibility of collapse of the simplex. See [2,4] for details.

Specifying Simplex Type

By default, a right simplex is initialized. If you prefer a regular simplex, you must call void ChooseRegularSimplex() before calling BeginSearch(). This causes a regular simplex to be initialized, using only the first edge length in the vector of edge lengths specified by the user (or an edge length of 2.0 by default). The method void ChooseRightSimplex() is also available. To read an arbitrary simplex in from a file, use the method void ReadInFile(istream& fp). The program main_arb (the source listing is main_arb.cc) shows an example of this. It optimizes Powell's quartic using the simplex in the input file arbitrary.dat. Be sure that arbitrary.dat is in the directory where you are running main_arb, or it won't be able to find it.

When preparing an input file, remember that C++ uses row-major order for doubly-indexed arrays. Our Matrix class uses such an array as its underlying data structure, so if you enter your coordinates by column rather than row, you will not have the simplex you were looking for. The simplex encoded in arbitrary.dat will give the following when printed out with PrintDesign(): 

Point: 2 2 0 1 
Value: 515
Point: -4 0 -2 -3 
Value: 287
Point: 0 0 0 7 
Value: 24255
Point: 0 0 0 -8 
Value: 41280
Point: 2 -5 -3 7 
Value: 9055

Specifying initial settings in the Nelder-Mead Search

The settings assigned by default for fields alpha, beta, and gamma are those recommended by Nelder and Mead [5], as specified by [3]. These were found by the authors to give the fastest convergence results compared to other strategies. The value of 0.5 for sigma was used throughout, and is the standard value for the shrinking coefficient in Simplex searches.

Depending on the problem, however, the user may wish to use different values for these variables. Most NMSearch constructors take these as parameters; examples of how to fill in the blanks can be found in listings for main_spec.cc and main_arb.cc. Alternatively, they can be changed one at a time using these mutator methods: 

   void SetAlpha(double newAlpha)
   void SetBeta(double newBeta)
   void SetGamma(double newGamma)
   void SetSigma(double newSigma)

This also applies to Hybrid_NMSearch, in which case the above values (except for sigma would apply only to the NMSearch phase of the search.

Appendix A: Settings for Fields in DirectSearch classes


Fields first declared in DirectSearch class

  1. Matrix $<$double$>$ * design

  2. long dimension

  3. Vector $<$double$>$ * minPoint

  4. long functionCalls

  5. long maxCalls

  6. double stoppingStepLength

  7. void* some_object

  8. void(*fcn_name)(long dim, Vector<double> &x,

    double & function, bool & success, void* an_object)

  9. bool exact_count

  10. int IDnumber


Fields first declared in PatternSearch class

  1. long patternLength

  2. double delta

  3. double initialStepLength


Fields first declared in SimplexSearch class

  1. double *simplexValues

  2. Vector$<$double$>$ *starting_edgeLengths

  3. double sigma

  4. long minIndex

  5. long replacementIndex

  6. Vector$<$double$>$ *centroid

  7. int toleranceHit

  8. bool SimplexSpecified

  9. bool Stop_on_std

  10. double delta

  11. Vector$<$double$>$ *scratch, *scratch2


Fields first declared in PatternSearch-derived concrete classes

HJSearch/EdHJSearch

  1. double step

  2. double factor

MultiCompassSearch

  1. long NumSearches

  2. IterationCap

  3. long stream

  4. Vector$<$double$>$ *random_start_pt

  5. Bounds bds


Fields first declared in SimplexSearch-derived concrete classes


NMSearch

  1. double alpha

  2. double beta

  3. double gamma

  4. double NMdelta

  5. long maxIndex

  6. Vector$<$double$>$ *reflectionPt

  7. double reflectionPtValue

  8. Vector$<$double$>$ *expansionPt

  9. double expansionPtValue

  10. Vector$<$double$>$ *contractionPt

  11. double contractionPtValue

  12. double maxPrimePtValue

  13. long maxPrimePtId

SHHSearch

  1. long *simplexAges

  2. Vector$<$double$>$ *reflectionPt

  3. double reflectionPtValue

SMDSearch

  1. long *simplexVBits

  2. long currentIndex

  3. long refCurrentIndex

  4. Matrix$<$double$>$ *refSimplex

  5. double *refSimplexValues

  6. long *refSimplexVBits

Appendix B: ID Numbers of the Searches.



Search type ID Number
DirectSearch 1000
PatternSearch 2000
CoordinateSearch 2100
CompassSearch 2200
MultiCompassSearch 2210
HJSearch 2300
EdHJSearch 2400
NLessSearch 2500
SimplexSearch 3000
SHHSearch 3100
NMSearch 3200
Hybrid_NMSearch.cc 3210
SMDSearch 3300

Acknowledgements

Some of the material in this manual was adapted from Chris Siefert's user's manual for MAPS (available at http://www.cs.wm.edu/$\tilde{ }$va/software/maps/manual.html) and Tony Padula's documentation for his krigifier software (available at http://www.cs.wm.edu/$\tilde{ }$va/software/krigifier/documentation). Adam Gurson and Liz Dolan wrote the original search software from which DirectSearch and its derived classes were developed. My advisor, Dr.Virginia Torczon, gave vast amounts of help and advice throughout the project. And, of course, we gratefully acknowledge the support of the National Science Foundation, which supported the research underlying this project under grant number CCR-9734044.

Bibliography

1
Elizabeth D. Dolan, Robert Michael Lewis, and Virginia J. Torczon.
On the local convergence properties of pattern search.
In review for release as an ICASE technical report. Submitted to SIAM Journal on Optimization.

2
Adam P. Gurson.
Simplex search behavior in nonlinear optimization.
Honors Thesis, Department of Computer Science, College of William & Mary, Williamsburg, Virginia 23187-8795, May 2000. Accepted with Highest Honors.

3
Jeffrey C. Lagarias, James A. Reeds, Margaret H. Wright, and Paul E. Wright.
Convergence properties of the Nelder-Mead simplex method in low dimensions.
SIAM Journal on Optimization, 9(1):112-147, 1998.

4
K.I.M. McKinnon.
Convergence of the Nelder-Mead simplex method to a nonstationary point.
SIAM Journal on Optimization, 9(1):148-158, December 1998.

5
J. A. Nelder and R. Mead.
A simplex method for function minimization.
The Computer Journal, 7(4):308-313, January 1965.

About this document ...

User's Manual for DirectSearch and its Derived Classes

This document was generated using the LaTeX2HTML translator Version 99.2beta6 (1.42)

Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
Copyright © 1997, 1998, 1999, Ross Moore, Mathematics Department, Macquarie University, Sydney.

The command line arguments were:
latex2html -antialias -no_navigation -white -split 0 -t DSmanual -image_type gif manual.tex

The translation was initiated by Anne Shepherd on 2001-05-18

Anne Shepherd 2001-05-18