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