src/work/akos/simann.h
changeset 1144 1cfabf245433
parent 1096 1cfb25ef14d2
child 1145 99c1aa395a58
equal deleted inserted replaced
7:37dcb90c4739 8:af93b547b958
     1 #ifndef LEMON_SIMANN_H
     1 #ifndef LEMON_SIMANN_H
     2 #define LEMON_SIMANN_H
     2 #define LEMON_SIMANN_H
       
     3 
       
     4 /// \ingroup experimental
       
     5 /// \file
       
     6 /// \brief Simulated annealing framework.
       
     7 /// \author Akos Ladanyi
     3 
     8 
     4 #include <cstdlib>
     9 #include <cstdlib>
     5 #include <cmath>
    10 #include <cmath>
     6 #include <lemon/time_measure.h>
    11 #include <lemon/time_measure.h>
     7 
    12 
     8 namespace lemon {
    13 namespace lemon {
     9 
    14 
       
    15 /// \addtogroup experimental
       
    16 /// @{
       
    17 
    10   const double INFTY = 1e24;
    18   const double INFTY = 1e24;
    11 
    19 
       
    20   /*! \brief Simulated annealing base class. */
    12   class SimAnnBase {
    21   class SimAnnBase {
    13   public:
    22   public:
    14     class Controller;
    23     class Controller;
    15   private:
    24   private:
       
    25     /*! Pointer to the controller. */
    16     Controller *controller;
    26     Controller *controller;
    17   protected:
    27   protected:
       
    28     /*! \brief Cost of the current solution. */
    18     double curr_cost;
    29     double curr_cost;
       
    30     /*! \brief Cost of the best solution. */
    19     double best_cost;
    31     double best_cost;
       
    32     /*! \brief Cost of the previous solution. */
    20     double prev_cost;
    33     double prev_cost;
       
    34     /*! \brief Cost of the solution preceding the previous one. */
    21     double prev_prev_cost;
    35     double prev_prev_cost;
    22 
    36 
    23     virtual void mutate() = 0;
    37     /*! \brief Step to a neighbouring state. */
    24     virtual void revert() = 0;
    38     virtual void mutate() {}
    25     virtual void saveAsBest() = 0;
    39     /*! \brief Reverts the last mutate(). */
    26   public:
    40     virtual void revert() {}
       
    41     /*! \brief Saves the current solution as the best one. */
       
    42     virtual void saveAsBest() {}
       
    43   public:
       
    44     /*! \brief Constructor. */
    27     SimAnnBase() {
    45     SimAnnBase() {
    28       best_cost = prev_cost = prev_prev_cost = INFTY;
    46       best_cost = prev_cost = prev_prev_cost = INFTY;
    29     }
    47     }
       
    48     /*! \brief Sets the controller class to use. */
    30     void setController(Controller &_controller) {
    49     void setController(Controller &_controller) {
    31       controller = &_controller;
    50       controller = &_controller;
    32       controller->setBase(this);
    51       controller->setBase(this);
    33     }
    52     }
       
    53     /*! \brief Returns the cost of the current solution. */
    34     double getCurrCost() const { return curr_cost; }
    54     double getCurrCost() const { return curr_cost; }
       
    55     /*! \brief Returns the cost of the previous solution. */
    35     double getPrevCost() const { return prev_cost; }
    56     double getPrevCost() const { return prev_cost; }
       
    57     /*! \brief Returns the cost of the best solution. */
    36     double getBestCost() const { return best_cost; }
    58     double getBestCost() const { return best_cost; }
       
    59     /*! \brief Starts the annealing process. */
    37     void run() {
    60     void run() {
    38       controller->init();
    61       controller->init();
    39       do {
    62       do {
    40         mutate();
    63         mutate();
    41         if (controller->accept()) {
    64         if (controller->accept()) {
    53     }
    76     }
    54 
    77 
    55     /*! \brief A base class for controllers. */
    78     /*! \brief A base class for controllers. */
    56     class Controller {
    79     class Controller {
    57     public:
    80     public:
       
    81       /*! \brief Pointer to the simulated annealing base class. */
    58       SimAnnBase *base;
    82       SimAnnBase *base;
       
    83       /*! \brief Initializes the controller. */
    59       virtual void init() {}
    84       virtual void init() {}
    60       /*! \brief This is called when a neighbouring state gets accepted. */
    85       /*! \brief This is called when a neighbouring state gets accepted. */
    61       virtual void acceptEvent() {}
    86       virtual void acceptEvent() {}
    62       /*! \brief This is called when the accepted neighbouring state's cost is
    87       /*! \brief This is called when the accepted neighbouring state's cost is
    63        *  less than the best found one's.
    88        *  less than the best found one's.
    64        */
    89        */
    65       virtual void improveEvent() {}
    90       virtual void improveEvent() {}
    66       /*! \brief This is called when a neighbouring state gets rejected. */
    91       /*! \brief This is called when a neighbouring state gets rejected. */
    67       virtual void rejectEvent() {}
    92       virtual void rejectEvent() {}
       
    93       /*! \brief Sets the simulated annealing base class to use. */
    68       virtual void setBase(SimAnnBase *_base) { base = _base; }
    94       virtual void setBase(SimAnnBase *_base) { base = _base; }
    69       /*! */
    95       /*! \brief Decides whether to continue the annealing process or not. */
    70       virtual bool next() = 0;
    96       virtual bool next() = 0;
    71       /*! */
    97       /*! \brief Decides whether to accept the current solution or not. */
    72       virtual bool accept() = 0;
    98       virtual bool accept() = 0;
    73     };
    99     };
    74   };
   100   };
    75 
   101 
       
   102   /*! \brief Simulated annealing class. */
    76   template <typename E>
   103   template <typename E>
    77   class SimAnn : public SimAnnBase {
   104   class SimAnn : public SimAnnBase {
    78   private:
   105   private:
       
   106     /*! \brief Pointer to the current entity. */
    79     E *curr_ent;
   107     E *curr_ent;
       
   108     /*! \brief Pointer to the best entity. */
    80     E *best_ent;
   109     E *best_ent;
    81   public:
   110   public:
       
   111     /*! \brief Constructor. */
    82     SimAnn() : SimAnnBase() {}
   112     SimAnn() : SimAnnBase() {}
       
   113     /*! \brief Sets the initial entity. */
    83     void setEntity(E &ent) {
   114     void setEntity(E &ent) {
    84       curr_ent = new E(ent);
   115       curr_ent = new E(ent);
    85       best_ent = new E(ent);
   116       best_ent = new E(ent);
    86       curr_cost = curr_ent->getCost();
   117       curr_cost = curr_ent->getCost();
    87     }
   118     }
       
   119     /*! \brief Returns the best found entity. */
    88     E getBestEntity() { return *best_ent; }
   120     E getBestEntity() { return *best_ent; }
       
   121     /*! \brief Step to a neighbouring state. */
    89     void mutate() {
   122     void mutate() {
    90       prev_prev_cost = prev_cost;
   123       prev_prev_cost = prev_cost;
    91       prev_cost = curr_cost;
   124       prev_cost = curr_cost;
    92       curr_ent->mutate();
   125       curr_ent->mutate();
    93       curr_cost = curr_ent->getCost();
   126       curr_cost = curr_ent->getCost();
    94     }
   127     }
       
   128     /*! \brief Reverts the last mutate(). */
    95     void revert() {
   129     void revert() {
    96       curr_ent->revert();
   130       curr_ent->revert();
    97       curr_cost = prev_cost;
   131       curr_cost = prev_cost;
    98       prev_cost = prev_prev_cost;
   132       prev_cost = prev_prev_cost;
    99     }
   133     }
       
   134     /*! \brief Saves the current solution as the best one. */
   100     void saveAsBest() {
   135     void saveAsBest() {
   101       delete(best_ent);
   136       delete(best_ent);
   102       best_ent = new E(*curr_ent);
   137       best_ent = new E(*curr_ent);
   103       best_cost = curr_cost;
   138       best_cost = curr_cost;
   104     }
   139     }
   105   };
   140   };
   106 
   141 
       
   142   /*! \brief Skeleton of an entity class. */
   107   class EntitySkeleton {
   143   class EntitySkeleton {
   108   public:
   144   public:
   109     /*! \return the cost of the entity */
   145     /*! \brief Returns the cost of the entity. */
   110     double getCost() { return 0.0; }
   146     double getCost() { return 0.0; }
   111     /*! \brief Makes a minor change to the entity. */
   147     /*! \brief Makes a minor change to the entity. */
   112     void mutate() {}
   148     void mutate() {}
   113     /*! \brief Restores the entity to its previous state i.e. reverts the
   149     /*! \brief Restores the entity to its previous state i.e. reverts the
   114      *  effects of the last mutate.
   150      *  effects of the last mutate().
   115      */
   151      */
   116     void revert() {}
   152     void revert() {}
   117   };
   153   };
   118 
   154 
   119   /*! \brief A simple controller for the simulated annealing class.
   155   /*! \brief A simple controller for the simulated annealing class. */
   120    *  \todo Find a way to set the various parameters.
       
   121    */
       
   122   class SimpleController : public SimAnnBase::Controller {
   156   class SimpleController : public SimAnnBase::Controller {
   123   public:
   157   public:
   124     long iter, last_impr, max_iter, max_no_impr;
   158     /*! \brief Number of iterations. */
   125     double temp, ann_fact;
   159     long iter;
   126     /*! \param _max_iter maximum number of iterations
   160     /*! \brief Number of iterations which did not improve the solution since
       
   161      *  the last improvement. */
       
   162     long last_impr;
       
   163     /*! \brief Maximum number of iterations. */
       
   164     long max_iter;
       
   165     /*! \brief Maximum number of iterations which do not improve the
       
   166      *  solution. */
       
   167     long max_no_impr;
       
   168     /*! \brief Temperature. */
       
   169     double temp;
       
   170     /*! \brief Annealing factor. */
       
   171     double ann_fact;
       
   172     /*! \brief Constructor.
       
   173      *  \param _max_iter maximum number of iterations
   127      *  \param _max_no_impr maximum number of consecutive iterations which do
   174      *  \param _max_no_impr maximum number of consecutive iterations which do
   128      *         not yield a better solution
   175      *         not yield a better solution
   129      *  \param _temp initial temperature
   176      *  \param _temp initial temperature
   130      *  \param _ann_fact annealing factor
   177      *  \param _ann_fact annealing factor
   131      */
   178      */
   134     max_iter(_max_iter), max_no_impr(_max_no_impr), temp(_temp),
   181     max_iter(_max_iter), max_no_impr(_max_no_impr), temp(_temp),
   135     ann_fact(_ann_fact) {}
   182     ann_fact(_ann_fact) {}
   136     void acceptEvent() {
   183     void acceptEvent() {
   137       iter++;
   184       iter++;
   138     }
   185     }
       
   186     /*! \brief This is called when the accepted neighbouring state's cost is
       
   187      *  less than the best found one's.
       
   188      */
   139     void improveEvent() {
   189     void improveEvent() {
   140       last_impr = iter;
   190       last_impr = iter;
   141     }
   191     }
       
   192     /*! \brief This is called when a neighbouring state gets rejected. */
   142     void rejectEvent() {
   193     void rejectEvent() {
   143       iter++;
   194       iter++;
   144     }
   195     }
       
   196     /*! \brief Decides whether to continue the annealing process or not. Also
       
   197      *  decreases the temperature. */
   145     bool next() {
   198     bool next() {
   146       temp *= ann_fact;
   199       temp *= ann_fact;
   147       bool quit = (iter > max_iter) || (iter - last_impr > max_no_impr);
   200       bool quit = (iter > max_iter) || (iter - last_impr > max_no_impr);
   148       return !quit;
   201       return !quit;
   149     }
   202     }
       
   203     /*! \brief Decides whether to accept the current solution or not. */
   150     bool accept() {
   204     bool accept() {
   151       double cost_diff = base->getPrevCost() - base->getCurrCost();
   205       double cost_diff = base->getPrevCost() - base->getCurrCost();
   152       if (cost_diff < 0.0) {
   206       if (cost_diff < 0.0) {
   153         bool ret = drand48() <= exp(cost_diff / temp);
   207         bool ret = drand48() <= exp(cost_diff / temp);
   154         return ret;
   208         return ret;
   159     }
   213     }
   160   };
   214   };
   161 
   215 
   162   /*! \brief A controller with preset running time for the simulated annealing
   216   /*! \brief A controller with preset running time for the simulated annealing
   163    *  class.
   217    *  class.
   164    *  \todo Find a better name.
       
   165    */
   218    */
   166   class AdvancedController : public SimAnnBase::Controller {
   219   class AdvancedController : public SimAnnBase::Controller {
   167   private:
   220   private:
   168     Timer timer;
   221     Timer timer;
   169     /*! \param time the elapsed time in seconds */
   222     /*! \param time the elapsed time in seconds */
   170     virtual double threshold(double time) {
   223     virtual double threshold(double time) {
   171       // 1 / log(x)
       
   172       /*
       
   173       static double xm = 5.0 / end_time;
       
   174       static double ym = start_threshold / (1 / log(1.2) - 1 / log(5.0 + 1.2));
       
   175       return ym * (1 / log(xm * time + 1.2) - 1 / log(5.0 + 1.2));
       
   176       */
       
   177       return (-1.0) * start_threshold / end_time * time + start_threshold;
   224       return (-1.0) * start_threshold / end_time * time + start_threshold;
   178     }
   225     }
   179   public:
   226   public:
   180     double alpha, beta, gamma;
   227     double alpha;
   181     double end_time, start_time;
   228     double beta;
       
   229     double gamma;
       
   230     double end_time;
       
   231     double start_time;
   182     double start_threshold;
   232     double start_threshold;
   183     double avg_cost;
   233     double avg_cost;
   184     double temp, ann_fact;
   234     double temp;
       
   235     double ann_fact;
   185     bool warmup;
   236     bool warmup;
   186     /*! \param _end_time running time in seconds
   237     /*! \brief Constructor.
       
   238      *  \param _end_time running time in seconds
   187      *  \param _alpha parameter used to calculate the running average
   239      *  \param _alpha parameter used to calculate the running average
   188      *  \param _beta parameter used to decrease the annealing factor
   240      *  \param _beta parameter used to decrease the annealing factor
   189      *  \param _gamma parameter used to increase the temperature
   241      *  \param _gamma parameter used to increase the temperature
   190      */
   242      */
   191     AdvancedController(double _end_time, double _alpha = 0.2,
   243     AdvancedController(double _end_time, double _alpha = 0.2,
   192     double _beta = 0.9, double _gamma = 1.6) : alpha(_alpha), beta(_beta),
   244     double _beta = 0.9, double _gamma = 1.6) : alpha(_alpha), beta(_beta),
   193     gamma(_gamma), end_time(_end_time), ann_fact(0.99999999), warmup(true) {}
   245     gamma(_gamma), end_time(_end_time), ann_fact(0.99999999), warmup(true) {}
   194     void init() {
   246     void init() {
   195       avg_cost = base->getCurrCost();
   247       avg_cost = base->getCurrCost();
   196     }
   248     }
       
   249     /*! \brief This is called when a neighbouring state gets accepted. */
   197     void acceptEvent() {
   250     void acceptEvent() {
   198       avg_cost = alpha * base->getCurrCost() + (1.0 - alpha) * avg_cost;
   251       avg_cost = alpha * base->getCurrCost() + (1.0 - alpha) * avg_cost;
   199       if (warmup) {
   252       if (warmup) {
   200         static int cnt = 0;
   253         static int cnt = 0;
   201         cnt++;
   254         cnt++;
   202         if (cnt >= 100) {
   255         if (cnt >= 100) {
   203           // calculate starting threshold and starting temperature
   256           // calculate starting threshold and starting temperature
   204           start_threshold = 5.0 * fabs(base->getBestCost() - avg_cost);
   257           start_threshold = 5.0 * fabs(base->getBestCost() - avg_cost);
   205           //temp = max_cost_diff / log(0.5);
       
   206           temp = 10000.0;
   258           temp = 10000.0;
   207           warmup = false;
   259           warmup = false;
   208           timer.reset();
   260           timer.reset();
   209         }
   261         }
   210       }
   262       }
   211     }
   263     }
   212     void improveEvent() {
   264     /*! \brief Decides whether to continue the annealing process or not. */
   213     }
       
   214     void rejectEvent() {
       
   215     }
       
   216     bool next() {
   265     bool next() {
   217       if (warmup) {
   266       if (warmup) {
   218         return true;
   267         return true;
   219       }
   268       }
   220       else {
   269       else {
   224           ann_fact *= beta;
   273           ann_fact *= beta;
   225         }
   274         }
   226         else {
   275         else {
   227           // increase the temperature
   276           // increase the temperature
   228           temp *= gamma;
   277           temp *= gamma;
   229           ann_fact = 0.99999999;
   278           ann_fact = 0.99999999; // !!!!!!!!!!!
   230         }
   279         }
   231         temp *= ann_fact;
   280         temp *= ann_fact;
   232         return elapsed_time < end_time;
   281         return elapsed_time < end_time;
   233       }
   282       }
   234     }
   283     }
       
   284     /*! \brief Decides whether to accept the current solution or not. */
   235     bool accept() {
   285     bool accept() {
   236       if (warmup) {
   286       if (warmup) {
   237         // we accept eveything during the "warm up" phase
   287         // we accept eveything during the "warm up" phase
   238         return true;
   288         return true;
   239       }
   289       }
   247         }
   297         }
   248       }
   298       }
   249     }
   299     }
   250   };
   300   };
   251 
   301 
       
   302 /// @}
       
   303 
   252 }
   304 }
   253 
   305 
   254 #endif
   306 #endif