smdsearch.h

在linux下面实现的单纯性算法的源代码

/*SMDSearch.h
 *Declarations of Sequential version of Torczon's Multi-Directional Search
 *Adam Gurson College of William & Mary 2000
 */

#ifndef _SMDSearch_
#define _SMDSearch_

#include <iostream.h>
#include <fstream.h>
#include <math.h>
#include <stdlib.h>
#include "objective.h"
#include "vec.h"
#include "cmat.h"
#include "subscrpt.h"


#define stoppingStepLength 10e-8 /*sqrt(fabs(3.0 * (4.0/3.0 - 1.0) - 1.0))*/
#ifndef DEBUG
#define DEBUG 0
#endif
//#ifndef maxCalls
//#define maxCalls 200
//#endif


class SMDSearch
{ 
 public:
 // constructors and destructors

   SMDSearch(int dim);
   // default constructor

   SMDSearch(int dim, double Sigma);
   // constructor which allows shrinking coefficient initialization

   SMDSearch(const SMDSearch& Original);
   // deep copy constructor

   ~SMDSearch();
   // destructor

 // algorithmic routines

   void ExploratoryMoves();
   // use SMD to find function minimum

   void ReplaceSimplexPoint(int index, const Vector<double>& newPoint);
   // replaces simplex point indexed at index with newPoint

   void CalculateFunctionValue(int index);
   // finds the f(x) value for the simplex point indexed at index and
   // replaces the proper value in simplexValues

   void SetSigma(double newSigma);
   // allows the user to set a new value for
   // the shrinking coefficient

   bool Stop();
   // returns true if the stopping criteria have been satisfied

   void fcnCall(int n, double *x, double& f, int& flag);
   // indirection of function call for purposes of keeping an accurate
   // tally of the number of function calls

 // Simplex-altering functions

   void InitRegularTriangularSimplex(const Vector<double> *basePoint,
                                     const double edgeLength);
   // deletes any existing simplex and replaces it with a regular
   // triangular simplex in the following manner:
   //
   // basePoint points to a point that will be the "origin" of the
   //    simplex points (it will be a part of the simplex and
   //    its function value is found here)
   // edgeLength is the length of each edge of the "triangle"
   //
   // functionCalls is reset to 0
   // delta is set to edgeLength
   //
   // NOTE: basePoint is assumed to be of proper dimension

   void InitFixedLengthRightSimplex(const Vector<double> *basePoint,
                                    const double edgeLength);
   // deletes any existing simplex and replaces it with a right-angle
   // simplex in the following manner:
   //
   // basePoint points to a point that will be the "origin" of the
   //    simplex points (it will be a part of the simplex and
   //    its function value is found here)
   // edgeLength is to be the length of each simplex side extending
   //    from the basePoint along each positive coordinate direction.
   //
   // functionCalls is reset to 0
   // delta is set to edgeLength
   //
   // NOTE: basePoint is assumed to be of proper dimension

   void InitVariableLengthRightSimplex(const Vector<double> *basePoint,
                                       const double* edgeLengths);
   // deletes any existing simplex and replaces it with a right-angle
   // simplex in the following manner:
   //
   // basePoint points to a point that will be the "origin" of the
   //    simplex points (it will be a part of the simplex and
   //    its function value is found here)
   // edgeLengths points to an array of n doubles, where n is the
   //    dimension of the given search. x_1 will then be located
   //    a distance of edgeLengths[0] away from the basepoint along the
   //    the x_1 axis, x_2 is edgeLengths[1] away on the x_2 axis, etc.
   //
   // functionCalls is reset to 0
   // delta is set to the largest value in edgeLength[]
   //
   // NOTE: basePoint and edgeLengths are assumed to be of proper dimension

   void InitGeneralSimplex(const Matrix<double> *plex);
   // deletes any existing simplex and replaces it with the one
   // pointed to by plex
   //
   // functionCalls is reset to 0
   // delta is set to the length of the longest simplex side
   //
   // NOTE: THIS ASSUMES THAT plex IS OF PROPER DIMENSION
   //       AND THAT basePoint is in the LAST row of plex
   //
   // The basePoint is a point that will be the "origin" of the
   //    simplex points (it will be a part of the simplex and
   //    its function value is calculated here)

   void ReadSimplexFile(istream& fp);
   // may also pass cin as input stream if desired
   // input the values of each trial point
   // NOTE: THIS FUNCTION WILL ONLY ACCEPT n+1 POINTS
   //
   // NOTE: assumes that the basePoint is the last point entered
   //
   // functionCalls is reset to 0
   // delta is set to the length of the longest simplex side

 // Query functions

   int GetFunctionCalls() const;
   // number of objective function evaluations

   void GetMinPoint(Vector<double>* &minimum) const;
   // simplex point which generates the best objective function
   // value found thus far
   // USER SHOULD PASS JUST A NULL POINTER, WITHOUT PREALLOCATED MEMORY

   double GetMinVal() const;
   // best objective function value found thus far

   void GetCurrentSimplex(Matrix<double>* &plex) const;
   // performs a deep copy of the simplex to a Matrix pointer
   // points to a newly allocated chunk of memory upon return
   // USER SHOULD PASS JUST A NULL POINTER, WITHOUT PREALLOCATED MEMORY

   void GetCurrentSimplexValues(double* &simValues) const;
   // performs a deep copy of the simplexValues array to a double pointer
   // points to a newly allocated chunk of memory upon return
   // USER SHOULD PASS JUST A NULL POINTER, WITHOUT PREALLOCATED MEMORY

   void GetCurrentSimplexVBits(int* &simVBits) const;
   // performs a deep copy of the simplexVBits array to a double pointer
   // points to a newly allocated chunk of memory upon return
   // USER SHOULD PASS JUST A NULL POINTER, WITHOUT PREALLOCATED MEMORY

   int GetVarNo() const;
   // returns the dimension of the problem

   int GetTolHit() const;
   // returns toleranceHit

   void printSimplex() const;
   // prints out the primary simplex points by row, 
   // their corresponding f(x) values, their validity status, 
   // and the number of function calls thus far

   void printRefSimplex() const;
   // prints out the reflection simplex points by row, 
   // their corresponding f(x) values, their validity status, 
   // and the number of function calls thus far

 private:

   void CreateRefSimplex();
   // creates a reflection simplex from the primary simplex

   void SwitchSimplices();
   // swaps the primary and reflection simplices

   void ShrinkSimplex();
   // this function goes through the primary simplex and reduces the
   // lengths of the edges adjacent to the best vertex

   int GetAnotherIndex(int& index, int*& validBits);
   // this is used to find another INVALID point in the simplex
   // if all points are VALID, returns 0, otherwise 1

   void CalculateRefFunctionValue(int index);
   // Like CalculateFunctionValue(), but for the Reflection Simplex
   // (A user should not directly manipulate the reflection simplex,
   // hence the private status of this function)


   int dimensions;                // the number of variables
                                  //    (the dimension of the problem)
   Matrix<double> *simplex;       // the current simplex
   double *simplexValues;         // their corresponding f(x) values
   int *simplexVBits;             // valid bits for the simplex values
   int currentIndex;              // used to step through the arrays

   Matrix<double> *refSimplex;       // the reflection simplex
   double *refSimplexValues;         // their corresponding f(x) values
   int *refSimplexVBits;             // valid bits for the simplex values
   int refCurrentIndex;              // used to step through the arrays

   Vector<double> *minPoint;         // the actual minimizer (i.e. best point seen)
   double minValue;               // f(minimizer)
   int minIndex;                  // the row index of minPoint in the main simplex
   double delta;                  // our stopping criterion
   double sigma;                  // shrinking coefficient
   long functionCalls;            // tally of number of function calls
   int toleranceHit;              // 1 if stop due to tolerance, 0 if funcCalls

   // the following vectors are simply extra storage space
   // that are used by functions that require a vector of
   // size = dimensions
   Vector<double> *scratch, *scratch2;
};


#endif