9 /// \todo Doc should be improved |
9 /// \todo Doc should be improved |
10 /// \author Akos Ladanyi |
10 /// \author Akos Ladanyi |
11 |
11 |
12 #include <cstdlib> |
12 #include <cstdlib> |
13 #include <cmath> |
13 #include <cmath> |
|
14 #include <limits> |
14 #include <lemon/time_measure.h> |
15 #include <lemon/time_measure.h> |
15 |
16 |
16 namespace lemon { |
17 namespace lemon { |
17 |
18 |
18 /// \addtogroup experimental |
19 /// \addtogroup experimental |
19 /// @{ |
20 /// @{ |
20 |
21 |
21 /*! \brief A base class for controllers. */ |
22 /// \brief A base class for controllers. |
22 class ControllerBase { |
23 class ControllerBase { |
|
24 public: |
23 friend class SimAnnBase; |
25 friend class SimAnnBase; |
24 public: |
26 /// \brief Pointer to the simulated annealing base class. |
25 /*! \brief Pointer to the simulated annealing base class. */ |
|
26 SimAnnBase *simann; |
27 SimAnnBase *simann; |
27 /*! \brief Initializes the controller. */ |
28 /// \brief Initializes the controller. |
28 virtual void init() {} |
29 virtual void init() {} |
29 /*! \brief This is called when a neighbouring state gets accepted. */ |
30 /// \brief This is called by the simulated annealing class when a |
|
31 /// neighbouring state gets accepted. |
30 virtual void acceptEvent() {} |
32 virtual void acceptEvent() {} |
31 /*! \brief This is called when the accepted neighbouring state's cost is |
33 /// \brief This is called by the simulated annealing class when the |
32 * less than the best found one's. |
34 /// accepted neighbouring state's cost is less than the best found one's. |
33 */ |
|
34 virtual void improveEvent() {} |
35 virtual void improveEvent() {} |
35 /*! \brief This is called when a neighbouring state gets rejected. */ |
36 /// \brief This is called by the simulated annealing class when a |
|
37 /// neighbouring state gets rejected. |
36 virtual void rejectEvent() {} |
38 virtual void rejectEvent() {} |
37 /*! \brief Decides whether to continue the annealing process or not. */ |
39 /// \brief Decides whether to continue the annealing process or not. |
38 virtual bool next() = 0; |
40 virtual bool next() = 0; |
39 /*! \brief Decides whether to accept the current solution or not. */ |
41 /// \brief Decides whether to accept the current solution or not. |
40 virtual bool accept() = 0; |
42 virtual bool accept() = 0; |
41 }; |
43 /// \brief Destructor. |
42 |
44 virtual ~ControllerBase() {} |
43 /*! \brief Skeleton of an entity class. */ |
45 }; |
|
46 |
|
47 /// \brief Skeleton of an entity class. |
44 class EntityBase { |
48 class EntityBase { |
45 public: |
49 public: |
46 /*! \brief Makes a minor change to the entity. |
50 /// \brief Makes a minor change to the entity. |
47 * \return the new cost |
51 /// \return the new cost |
48 */ |
|
49 virtual double mutate() = 0; |
52 virtual double mutate() = 0; |
50 /*! \brief Restores the entity to its previous state i.e. reverts the |
53 /// \brief Restores the entity to its previous state i.e. reverts the |
51 * effects of the last mutate(). |
54 /// effects of the last mutate(). |
52 */ |
|
53 virtual void revert() = 0; |
55 virtual void revert() = 0; |
54 /*! \brief Makes a copy of the entity. */ |
56 /// \brief Makes a copy of the entity. |
55 virtual EntityBase* clone() = 0; |
57 virtual EntityBase* clone() = 0; |
56 /*! \brief Makes a major change to the entity. */ |
58 /// \brief Makes a major change to the entity. |
57 virtual void randomize() = 0; |
59 virtual void randomize() = 0; |
58 }; |
60 /// \brief Destructor. |
59 |
61 virtual ~EntityBase() {} |
60 /*! \brief Simulated annealing base class. */ |
62 }; |
|
63 |
|
64 /// \brief Simulated annealing abstract base class. |
|
65 /// Can be used to derive a custom simulated annealing class if \ref SimAnn |
|
66 /// doesn't fit your needs. |
61 class SimAnnBase { |
67 class SimAnnBase { |
62 private: |
68 private: |
63 /*! Pointer to the controller. */ |
69 /// \brief Pointer to the controller. |
64 ControllerBase *controller; |
70 ControllerBase *controller; |
65 /*! \brief Cost of the current solution. */ |
71 /// \brief Cost of the current solution. |
66 double curr_cost; |
72 double curr_cost; |
67 /*! \brief Cost of the best solution. */ |
73 /// \brief Cost of the best solution. |
68 double best_cost; |
74 double best_cost; |
69 /*! \brief Cost of the previous solution. */ |
75 /// \brief Cost of the previous solution. |
70 double prev_cost; |
76 double prev_cost; |
71 /*! \brief Cost of the solution preceding the previous one. */ |
77 /// \brief Cost of the solution preceding the previous one. |
72 double prev_prev_cost; |
78 double prev_prev_cost; |
73 /*! \brief Number of iterations. */ |
79 /// \brief Number of iterations. |
74 long iter; |
80 long iter; |
75 /*! \brief Number of iterations which did not improve the solution since |
81 /// \brief Number of iterations which did not improve the solution since |
76 * the last improvement. */ |
82 /// the last improvement. |
77 long last_impr; |
83 long last_impr; |
78 protected: |
84 protected: |
79 /*! \brief Step to a neighbouring state. */ |
85 /// \brief Step to a neighbouring state. |
80 virtual double mutate() = 0; |
86 virtual double mutate() = 0; |
81 /*! \brief Reverts the last mutate(). */ |
87 /// \brief Reverts the last mutate(). |
82 virtual void revert() = 0; |
88 virtual void revert() = 0; |
83 /*! \brief Saves the current solution as the best one. */ |
89 /// \brief Saves the current solution as the best one. |
84 virtual void saveAsBest() = 0; |
90 virtual void saveAsBest() = 0; |
85 /*! \brief Does initializations before each run. */ |
91 /// \brief Does initializations before each run. |
86 virtual void init() { |
92 virtual void init() { |
87 controller->init(); |
93 controller->init(); |
88 curr_cost = prev_cost = prev_prev_cost = best_cost = |
94 curr_cost = prev_cost = prev_prev_cost = best_cost = |
89 std::numeric_limits<double>::infinity(); |
95 std::numeric_limits<double>::infinity(); |
90 iter = last_impr = 0; |
96 iter = last_impr = 0; |
91 } |
97 } |
92 public: |
98 public: |
93 /*! \brief Sets the controller class to use. */ |
99 /// \brief Sets the controller class to use. |
94 void setController(ControllerBase &_controller) { |
100 void setController(ControllerBase &_controller) { |
95 controller = &_controller; |
101 controller = &_controller; |
96 controller->simann = this; |
102 controller->simann = this; |
97 } |
103 } |
98 /*! \brief Returns the cost of the current solution. */ |
104 /// \brief Returns the cost of the current solution. |
99 double getCurrCost() const { return curr_cost; } |
105 double getCurrCost() const { return curr_cost; } |
100 /*! \brief Returns the cost of the previous solution. */ |
106 /// \brief Returns the cost of the previous solution. |
101 double getPrevCost() const { return prev_cost; } |
107 double getPrevCost() const { return prev_cost; } |
102 /*! \brief Returns the cost of the best solution. */ |
108 /// \brief Returns the cost of the best solution. |
103 double getBestCost() const { return best_cost; } |
109 double getBestCost() const { return best_cost; } |
104 /*! \brief Returns the number of iterations. */ |
110 /// \brief Returns the number of iterations done. |
105 long getIter() const { return iter; } |
111 long getIter() const { return iter; } |
106 /*! \brief Returns the number of the last iteration when the solution was |
112 /// \brief Returns the ordinal number of the last iteration when the |
107 * improved. |
113 /// solution was improved. |
108 */ |
|
109 long getLastImpr() const { return last_impr; } |
114 long getLastImpr() const { return last_impr; } |
110 /*! \brief Performs one iteration. */ |
115 /// \brief Performs one iteration. |
111 bool step() { |
116 bool step() { |
112 iter++; |
117 iter++; |
113 prev_prev_cost = prev_cost; |
118 prev_prev_cost = prev_cost; |
114 prev_cost = curr_cost; |
119 prev_cost = curr_cost; |
115 curr_cost = mutate(); |
120 curr_cost = mutate(); |
128 prev_cost = prev_prev_cost; |
133 prev_cost = prev_prev_cost; |
129 controller->rejectEvent(); |
134 controller->rejectEvent(); |
130 } |
135 } |
131 return controller->next(); |
136 return controller->next(); |
132 } |
137 } |
133 /*! \brief Performs a given number of iterations. |
138 /// \brief Performs a given number of iterations. |
134 * \param n the number of iterations |
139 /// \param n the number of iterations |
135 */ |
|
136 bool step(int n) { |
140 bool step(int n) { |
137 for(; n > 0 && step(); --n) ; |
141 for(; n > 0 && step(); --n) ; |
138 return !n; |
142 return !n; |
139 } |
143 } |
140 /*! \brief Starts the annealing process. */ |
144 /// \brief Starts the annealing process. |
141 void run() { |
145 void run() { |
142 init(); |
146 init(); |
143 do { } while (step()); |
147 do { } while (step()); |
144 } |
148 } |
145 }; |
149 /// \brief Destructor. |
146 |
150 virtual ~SimAnnBase() {} |
147 /*! \brief Simulated annealing class. */ |
151 }; |
|
152 |
|
153 /// \brief Simulated annealing class. |
148 class SimAnn : public SimAnnBase { |
154 class SimAnn : public SimAnnBase { |
149 private: |
155 private: |
150 /*! \brief Pointer to the current entity. */ |
156 /// \brief Pointer to the current entity. |
151 EntityBase *curr_ent; |
157 EntityBase *curr_ent; |
152 /*! \brief Pointer to the best entity. */ |
158 /// \brief Pointer to the best entity. |
153 EntityBase *best_ent; |
159 EntityBase *best_ent; |
154 /*! \brief Does initializations before each run. */ |
160 /// \brief Does initializations before each run. |
155 void init() { |
161 void init() { |
156 SimAnnBase::init(); |
162 SimAnnBase::init(); |
157 if (best_ent) delete best_ent; |
163 if (best_ent) delete best_ent; |
158 best_ent = NULL; |
164 best_ent = NULL; |
159 curr_ent->randomize(); |
165 curr_ent->randomize(); |
160 } |
166 } |
161 public: |
167 public: |
162 /*! \brief Constructor. */ |
168 /// \brief Constructor. |
163 SimAnn() : curr_ent(NULL), best_ent(NULL) {} |
169 SimAnn() : curr_ent(NULL), best_ent(NULL) {} |
164 /*! \brief Destructor. */ |
170 /// \brief Destructor. |
165 virtual ~SimAnn() { |
171 virtual ~SimAnn() { |
166 if (best_ent) delete best_ent; |
172 if (best_ent) delete best_ent; |
167 } |
173 } |
168 /*! \brief Step to a neighbouring state. */ |
174 /// \brief Step to a neighbouring state. |
169 double mutate() { |
175 double mutate() { |
170 return curr_ent->mutate(); |
176 return curr_ent->mutate(); |
171 } |
177 } |
172 /*! \brief Reverts the last mutate(). */ |
178 /// \brief Reverts the last mutate(). |
173 void revert() { |
179 void revert() { |
174 curr_ent->revert(); |
180 curr_ent->revert(); |
175 } |
181 } |
176 /*! \brief Saves the current solution as the best one. */ |
182 /// \brief Saves the current solution as the best one. |
177 void saveAsBest() { |
183 void saveAsBest() { |
178 if (best_ent) delete best_ent; |
184 if (best_ent) delete best_ent; |
179 best_ent = curr_ent->clone(); |
185 best_ent = curr_ent->clone(); |
180 } |
186 } |
181 /*! \brief Sets the current entity. */ |
187 /// \brief Sets the current entity. |
182 void setEntity(EntityBase &_ent) { |
188 void setEntity(EntityBase &_ent) { |
183 curr_ent = &_ent; |
189 curr_ent = &_ent; |
184 } |
190 } |
185 /*! \brief Returns a copy of the best found entity. */ |
191 /// \brief Returns a copy of the best found entity. |
186 EntityBase* getBestEntity() { return best_ent->clone(); } |
192 EntityBase* getBestEntity() { return best_ent->clone(); } |
187 }; |
193 }; |
188 |
194 |
189 /*! \brief A simple controller for the simulated annealing class. */ |
195 /// \brief A simple controller for the simulated annealing class. |
|
196 /// This controller starts from a given initial temperature and evenly |
|
197 /// decreases it. |
190 class SimpleController : public ControllerBase { |
198 class SimpleController : public ControllerBase { |
191 public: |
199 private: |
192 /*! \brief Maximum number of iterations. */ |
200 /// \brief Maximum number of iterations. |
193 long max_iter; |
201 long max_iter; |
194 /*! \brief Maximum number of iterations which do not improve the |
202 /// \brief Maximum number of iterations which do not improve the |
195 * solution. */ |
203 /// solution. |
196 long max_no_impr; |
204 long max_no_impr; |
197 /*! \brief Temperature. */ |
205 /// \brief Temperature. |
198 double temp; |
206 double temp; |
199 /*! \brief Annealing factor. */ |
207 /// \brief Annealing factor. |
200 double ann_fact; |
208 double ann_fact; |
201 /*! \brief Constructor. |
209 /// \brief Constructor. |
202 * \param _max_iter maximum number of iterations |
210 /// \param _max_iter maximum number of iterations |
203 * \param _max_no_impr maximum number of consecutive iterations which do |
211 /// \param _max_no_impr maximum number of consecutive iterations which do |
204 * not yield a better solution |
212 /// not yield a better solution |
205 * \param _temp initial temperature |
213 /// \param _temp initial temperature |
206 * \param _ann_fact annealing factor |
214 /// \param _ann_fact annealing factor |
207 */ |
215 public: |
208 SimpleController(long _max_iter = 500000, long _max_no_impr = 20000, |
216 SimpleController(long _max_iter = 500000, long _max_no_impr = 20000, |
209 double _temp = 1000.0, double _ann_fact = 0.9999) : max_iter(_max_iter), |
217 double _temp = 1000.0, double _ann_fact = 0.9999) : max_iter(_max_iter), |
210 max_no_impr(_max_no_impr), temp(_temp), ann_fact(_ann_fact) |
218 max_no_impr(_max_no_impr), temp(_temp), ann_fact(_ann_fact) |
211 { |
219 { |
212 srand48(time(0)); |
220 srand48(time(0)); |
213 } |
221 } |
214 /*! \brief This is called when a neighbouring state gets accepted. */ |
222 /// \brief This is called when a neighbouring state gets accepted. |
215 void acceptEvent() {} |
223 void acceptEvent() {} |
216 /*! \brief This is called when the accepted neighbouring state's cost is |
224 /// \brief This is called when the accepted neighbouring state's cost is |
217 * less than the best found one's. |
225 /// less than the best found one's. |
218 */ |
|
219 void improveEvent() {} |
226 void improveEvent() {} |
220 /*! \brief This is called when a neighbouring state gets rejected. */ |
227 /// \brief This is called when a neighbouring state gets rejected. |
221 void rejectEvent() {} |
228 void rejectEvent() {} |
222 /*! \brief Decides whether to continue the annealing process or not. Also |
229 /// \brief Decides whether to continue the annealing process or not. Also |
223 * decreases the temperature. */ |
230 /// decreases the temperature. |
224 bool next() { |
231 bool next() { |
225 temp *= ann_fact; |
232 temp *= ann_fact; |
226 bool quit = (simann->getIter() > max_iter) || |
233 bool quit = (simann->getIter() > max_iter) || |
227 (simann->getIter() - simann->getLastImpr() > max_no_impr); |
234 (simann->getIter() - simann->getLastImpr() > max_no_impr); |
228 return !quit; |
235 return !quit; |
229 } |
236 } |
230 /*! \brief Decides whether to accept the current solution or not. */ |
237 /// \brief Decides whether to accept the current solution or not. |
231 bool accept() { |
238 bool accept() { |
232 double cost_diff = simann->getPrevCost() - simann->getCurrCost(); |
239 double cost_diff = simann->getCurrCost() - simann->getPrevCost(); |
233 return (drand48() <= exp(cost_diff / temp)); |
240 return (drand48() <= exp(-(cost_diff / temp))); |
234 } |
241 } |
235 }; |
242 /// \brief Destructor. |
236 |
243 virtual ~SimpleController() {} |
237 /*! \brief A controller with preset running time for the simulated annealing |
244 }; |
238 * class. |
245 |
239 * |
246 /// \brief A controller with preset running time for the simulated annealing |
240 * With this controller you can set the running time of the annealing |
247 /// class. |
241 * process in advance. It works the following way: the controller measures |
248 /// With this controller you can set the running time of the annealing |
242 * a kind of divergence. The divergence is the difference of the average |
249 /// process in advance. It works the following way: the controller measures |
243 * cost of the recently found solutions the cost of the best found one. In |
250 /// a kind of divergence. The divergence is the difference of the average |
244 * case this divergence is greater than a given threshold, then we decrease |
251 /// cost of the recently found solutions the cost of the best found one. In |
245 * the annealing factor, that is we cool the system faster. In case the |
252 /// case this divergence is greater than a given threshold, then we decrease |
246 * divergence is lower than the threshold, then we increase the temperature. |
253 /// the annealing factor, that is we cool the system faster. In case the |
247 * The threshold is a function of the elapsed time which reaches zero at the |
254 /// divergence is lower than the threshold, then we increase the temperature. |
248 * desired end time. |
255 /// The threshold is a function of the elapsed time which reaches zero at the |
249 */ |
256 /// desired end time. |
250 class AdvancedController : public ControllerBase { |
257 class AdvancedController : public ControllerBase { |
251 private: |
258 private: |
|
259 /// \brief Timer class to measure the elapsed time. |
252 Timer timer; |
260 Timer timer; |
253 /*! \param time the elapsed time in seconds */ |
261 /// \brief Calculates the threshold value. |
|
262 /// \param time the elapsed time in seconds |
254 virtual double threshold(double time) { |
263 virtual double threshold(double time) { |
255 return (-1.0) * start_threshold / end_time * time + start_threshold; |
264 return (-1.0) * start_threshold / end_time * time + start_threshold; |
256 } |
265 } |
257 public: |
266 /// \brief Parameter used to calculate the running average. |
258 double alpha; |
267 double alpha; |
|
268 /// \brief Parameter used to decrease the annealing factor. |
259 double beta; |
269 double beta; |
|
270 /// \brief Parameter used to increase the temperature. |
260 double gamma; |
271 double gamma; |
261 /*! \brief The time at the end of the algorithm. */ |
272 /// \brief The time at the end of the algorithm. |
262 double end_time; |
273 double end_time; |
263 /*! \brief The time at the start of the algorithm. */ |
274 /// \brief The time at the start of the algorithm. |
264 double start_time; |
275 double start_time; |
265 /*! \brief Starting threshold. */ |
276 /// \brief Starting threshold. |
266 double start_threshold; |
277 double start_threshold; |
267 /*! \brief Average cost of recent solutions. */ |
278 /// \brief Average cost of recent solutions. |
268 double avg_cost; |
279 double avg_cost; |
269 /*! \brief Temperature. */ |
280 /// \brief Temperature. |
270 double temp; |
281 double temp; |
271 /*! \brief Annealing factor. */ |
282 /// \brief Annealing factor. |
272 double ann_fact; |
283 double ann_fact; |
273 /*! \brief Initial annealing factor. */ |
284 /// \brief Initial annealing factor. |
274 double init_ann_fact; |
285 double init_ann_fact; |
275 bool warmup; |
286 /// \brief True when the annealing process has been started. |
276 /*! \brief Constructor. |
287 bool start; |
277 * \param _end_time running time in seconds |
288 public: |
278 * \param _alpha parameter used to calculate the running average |
289 /// \brief Constructor. |
279 * \param _beta parameter used to decrease the annealing factor |
290 /// \param _end_time running time in seconds |
280 * \param _gamma parameter used to increase the temperature |
291 /// \param _alpha parameter used to calculate the running average |
281 * \param _ann_fact initial annealing factor |
292 /// \param _beta parameter used to decrease the annealing factor |
282 */ |
293 /// \param _gamma parameter used to increase the temperature |
|
294 /// \param _ann_fact initial annealing factor |
283 AdvancedController(double _end_time, double _alpha = 0.2, |
295 AdvancedController(double _end_time, double _alpha = 0.2, |
284 double _beta = 0.9, double _gamma = 1.6, double _ann_fact = 0.9999) : |
296 double _beta = 0.9, double _gamma = 1.6, double _ann_fact = 0.9999) : |
285 alpha(_alpha), beta(_beta), gamma(_gamma), end_time(_end_time), |
297 alpha(_alpha), beta(_beta), gamma(_gamma), end_time(_end_time), |
286 ann_fact(_ann_fact), init_ann_fact(_ann_fact), warmup(true) |
298 ann_fact(_ann_fact), init_ann_fact(_ann_fact), start(false) |
287 { |
299 { |
288 srand48(time(0)); |
300 srand48(time(0)); |
289 } |
301 } |
|
302 /// \brief Does initializations before each run. |
290 void init() { |
303 void init() { |
291 avg_cost = simann->getCurrCost(); |
304 avg_cost = simann->getCurrCost(); |
292 } |
305 } |
293 /*! \brief This is called when a neighbouring state gets accepted. */ |
306 /// \brief This is called when a neighbouring state gets accepted. |
294 void acceptEvent() { |
307 void acceptEvent() { |
295 avg_cost = alpha * simann->getCurrCost() + (1.0 - alpha) * avg_cost; |
308 avg_cost = alpha * simann->getCurrCost() + (1.0 - alpha) * avg_cost; |
296 if (warmup) { |
309 if (!start) { |
297 static int cnt = 0; |
310 static int cnt = 0; |
298 cnt++; |
311 cnt++; |
299 if (cnt >= 100) { |
312 if (cnt >= 100) { |
300 // calculate starting threshold and starting temperature |
313 // calculate starting threshold and starting temperature |
301 start_threshold = 5.0 * fabs(simann->getBestCost() - avg_cost); |
314 start_threshold = 5.0 * fabs(simann->getBestCost() - avg_cost); |
302 temp = 10000.0; |
315 temp = 10000.0; |
303 warmup = false; |
316 start = true; |
304 timer.restart(); |
317 timer.restart(); |
305 } |
318 } |
306 } |
319 } |
307 } |
320 } |
308 /*! \brief Decides whether to continue the annealing process or not. */ |
321 /// \brief Decides whether to continue the annealing process or not. |
309 bool next() { |
322 bool next() { |
310 if (warmup) { |
323 if (!start) { |
311 return true; |
324 return true; |
312 } |
325 } |
313 else { |
326 else { |
314 double elapsed_time = timer.getRealTime(); |
327 double elapsed_time = timer.realTime(); |
315 if (fabs(avg_cost - simann->getBestCost()) > threshold(elapsed_time)) { |
328 if (fabs(avg_cost - simann->getBestCost()) > threshold(elapsed_time)) { |
316 // decrease the annealing factor |
329 // decrease the annealing factor |
317 ann_fact *= beta; |
330 ann_fact *= beta; |
318 } |
331 } |
319 else { |
332 else { |