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