Main Page   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Namespace Members   Compound Members   File Members  

Optimizer.h

Go to the documentation of this file.
00001 /* WARANTY NOTICE AND COPYRIGHT
00002 This program is free software; you can redistribute it and/or
00003 modify it under the terms of the GNU General Public License
00004 as published by the Free Software Foundation; either version 2
00005 of the License, or (at your option) any later version.
00006 
00007 This program is distributed in the hope that it will be useful,
00008 but WITHOUT ANY WARRANTY; without even the implied warranty of
00009 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00010 GNU General Public License for more details.
00011 
00012 You should have received a copy of the GNU General Public License
00013 along with this program; if not, write to the Free Software
00014 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
00015 
00016 Copyright (C) Michael J. Meyer
00017 
00018 matmjm@mindspring.com
00019 spyqqqdia@yahoo.com
00020 
00021 */
00022 
00023 
00024 
00025 #ifndef martingale_optimizer_h    
00026 #define martingale_optimizer_h
00027 
00028 #include "TypedefsMacros.h"
00029 #include "Array.h"                            // direct members
00030 
00031 MTGL_BEGIN_NAMESPACE(Martingale)
00032 
00033 
00034 // we are using
00035 class LowDiscrepancySequence;
00036 
00037 
00038 
00059 /*******************************************************************************
00060  *
00061  *                     OPTIMIZER INTERFACE
00062  *
00063  ******************************************************************************/
00064  
00065 
00066 
00067     
00070 class Optimizer {
00071         
00072 protected:
00073     
00074     int n;                 // number of variables
00075 
00076 
00077 public:
00078                 
00080         virtual Real f(const RealArray1D& x) = 0;
00081     
00082     
00087     virtual const RealArray1D& search() = 0;
00088     
00089     
00092     Optimizer(int n_) : n(n_) { }
00093         virtual ~Optimizer(){ }
00094         
00095     
00096 }; // end Optimizer
00097 
00098 
00099 
00100    
00101 /*******************************************************************************
00102  *
00103  *                              DOWNHILL SIMPLEX
00104  *
00105  ******************************************************************************/
00106 
00107 
00112 class DownhillSimplex : public Optimizer {
00113     
00114   
00115    Array1D<RealArray1D*> vertices;         // vertices[i], i=0,...,n, are the 
00116                                            // vertices of the simplex
00117    RealArray1D sum;              // sum of all vertices
00118    RealArray1D newVertex;        // next vertex
00119    RealArray1D y;                // function values at the vertices
00120    
00121    Real fx;      // function at current point
00122    
00123    int min,      // index of minimum vertex
00124        max,      // index of vertex with highest function value
00125        max2,     // index of second highest vertex
00126        nStep;    // number of steps
00127    
00128    // The barycenter is updated dynamically whenever a vertex changes
00129    // since this is cheaper than doing a full computation in each step.
00130    
00131    bool verbose;
00132     
00133 public:
00134 
00135    RealArray1D& vertex(int i){ return *(vertices[i]); }
00136    Real& vertex(int i, int j){ return (*(vertices[i]))[j]; }
00137 
00138 // CONSTRUCTOR
00139     
00145     DownhillSimplex(const RealArray1D& x, Real delta, int steps, bool vbose);
00146    
00147     ~DownhillSimplex();
00148         
00149         
00154         void setInitialConditions();
00155 
00160     const RealArray1D& search();
00161 
00162 
00163          
00164 private:
00165  
00166          
00167 // REFLECTING, CONTRACTING, KEEPING TRACK OF HIGH-LOW POINTS     
00168     
00172     void setHiLowIndices();
00173         
00174         
00190     Real reflectVertex(int i, Real k);
00191         
00192         
00199     void contract(int i);
00200         
00201     
00209     void replaceWorst(Real w);
00210         
00211             
00212 }; // end DownhillSimplex
00213 
00214 
00215         
00219 class ConcreteDownhillSimplex : public DownhillSimplex {
00220         
00221          Real (*of)(const RealArray1D& x);   // pointer to the objective function
00222         
00223 public:
00224         
00227          ConcreteDownhillSimplex
00228          (Real (*f)(const RealArray1D&), RealArray1D& x, Real delta, int steps, bool verbose) :
00229      DownhillSimplex(x,delta,steps,verbose),
00230          of(f)
00231      {  
00232                  // now we can call f
00233                  setInitialConditions();
00234      }
00235          
00236          Real f(const RealArray1D& x){ return (*of)(x); }
00237 
00238 }; 
00239 
00240 
00241 
00242 /*******************************************************************************
00243  *
00244  *                        BFGS
00245  *
00246  ******************************************************************************/
00247 
00248 
00287 class BFGS : public Optimizer {
00288         
00289 protected:
00290         
00291     // some constants 
00292     static Real
00294     EPSX,
00296     KF,
00298     EPSG,
00300     EPS;
00301 
00302    
00303     Real stepmax;                 // maximum steplength in the line search
00304     int nVals,                    // maximum number of function evaluations
00305         fVals,                    // current number of function evaluations
00306             restarts,                 // number of times the optimizer has been restarted
00307             nRestarts;                // number of restarts budgeted
00308 
00309     // keep current
00310     RealArray1D 
00311         
00312                  x0,                  // the last point
00313              grad0,               // gradient at x0
00314              x,                   // the current point
00315              xDelta,              // x-x0
00316              grad,                // gradient at the current point x
00317              z,                   // next point to be computed (workspace)
00318              d,                   // current line direction
00319              gDelta,              // gradient difference grad-grad0
00320              hdg,                 // H(grad-grad0)
00321              u,                   // the additional vector in the bfgs update
00322              h;                   // directional increments to compute gradient
00323                          
00324     Real     fx0,                 // f(x0), function at the last point
00325              fx;                  // f(x), function at the current point
00326 
00327                  
00328     RealArray2D H;              // approximate inverse Hessian at the
00329                                   // current point.
00330     
00331     bool verbose;
00332 
00333     
00334 public:
00335 
00336 // ACCESSORS
00337 
00339     int getDimension() const { return n; }
00340     
00342     const RealArray1D& getGrad() const { return grad; }
00343     
00345     const RealArray1D& getX() const { return x; }
00346     
00348     const RealArray1D& getD() const { return d; }
00349         
00350 
00351 // CONSTRUCTOR
00352         
00353 
00354 
00363     BFGS(const RealArray1D& u, int vals, Real maxstep, 
00364              const RealArray1D& k, int resets=2, bool vbose=false);
00365         
00366                         
00371     void setInitialConditions();
00372     
00373                 
00374            
00375 // MIN SEARCH                
00376       
00382      const RealArray1D& search();
00383 
00384         
00385         
00386     
00387 private:
00388                 
00389 
00390         
00395     Real relativeSize(RealArray1D& x, RealArray1D& y);
00396         
00397     
00401     bool gradientIsZero();
00402     
00403     
00406     void bfgsUpdate();
00407         
00408         
00411         void resetHessian();
00412         
00413         
00414 // GRADIENT COMPUTATION  
00415     
00432     const RealArray1D& gradF(RealArray1D& x, Real fx);
00433         
00434    
00450     const RealArray1D& gradcdF(RealArray1D& x);
00451         
00452         
00453 // BACK TRACKING DURING LINE SEARCH                 
00454     
00465     Real backTrack(Real t1, Real t2, Real f1, Real f2, Real k);
00466 
00467 
00468 //  LINE SEARCH
00469   
00470     
00483     void lineSearch();
00484 
00485         
00486             
00487 }; // end BFGS
00488 
00489 
00490         
00494 class ConcreteBFGS : public BFGS {
00495         
00496          Real (*of)(const RealArray1D& x);   // pointer to the objective function
00497         
00498 public:
00499         
00502          ConcreteBFGS
00503          (Real (*f)(const RealArray1D&), RealArray1D& x, int nVals, 
00504       Real stepmax, RealArray1D& h, int nRestarts=3, bool verbose=false) :
00505      BFGS(x,nVals,stepmax,h,nRestarts,verbose),
00506          of(f)
00507      {   
00508                  // now we can call f
00509                  setInitialConditions();
00510          }
00511          
00512          Real f(const RealArray1D& x){ return (*of)(x); }
00513 
00514 }; 
00515 
00516 
00517 
00518 
00519 /*******************************************************************************
00520  *
00521  *              SOBOL SEARCH 
00522  *
00523  ******************************************************************************/
00524 
00532 class SobolSearch: public Optimizer {
00533         
00534 protected:
00535         
00536         int nPoints;                       // total number of search points
00537         Real q;                            // window contraction factor
00538         RealArray1D xOpt;                  // current best point
00539         RealArray1D x;                     // current point
00540         RealArray1D d;                     // window: all u with x_j-d_j < u_j < x_j+d_j
00541         
00542         bool verbose;
00543         
00544         LowDiscrepancySequence* lds;
00545         
00546         
00547 public:
00548         
00554         SobolSearch(const RealArray1D& x0, int nVals, const RealArray1D& delta, bool vbose=false);
00555         
00559         virtual bool isInDomain(const RealArray1D& u) const { return true; }
00560         
00561         const RealArray1D& search();
00562         
00563 }; // end SobolSearch
00564         
00565         
00566         
00567         
00571 class ConcreteSobolSearch : public SobolSearch {
00572         
00573          Real (*of)(const RealArray1D&);   // pointer to the objective function
00574         
00575 public:
00576         
00579          ConcreteSobolSearch
00580          (Real (*f)(const RealArray1D&), const RealArray1D& x, 
00581           int nVals, const RealArray1D& delta, bool verbose=false) :
00582      SobolSearch(x,nVals,delta,verbose),
00583          of(f)
00584      {    }
00585          
00586          Real f(const RealArray1D& x){ return (*of)(x); }
00587 
00588 };     
00589                 
00590         
00591         
00592         
00593 // TEST FUNCTIONS
00594 
00596 namespace ObjectiveFunction {
00597         
00598         Real function_1(const RealArray1D&);
00599         
00600 } // end namespace
00601 
00602 
00603 
00604 MTGL_END_NAMESPACE(Martingale)
00605 
00606 #endif

Generated on Mon Sep 22 02:16:32 2003 for Libor-Library by doxygen1.3-rc3