COIN-OR::LEMON - Graph Library

source: lemon-0.x/src/work/akos/simann.h @ 1139:f59038affc7e

Last change on this file since 1139:f59038affc7e was 1096:1cfb25ef14d2, checked in by Akos Ladanyi, 19 years ago

Various changes.

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