|
1 /* -*- C++ -*- |
|
2 * |
|
3 * This file is a part of LEMON, a generic C++ optimization library |
|
4 * |
|
5 * Copyright (C) 2003-2008 |
|
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport |
|
7 * (Egervary Research Group on Combinatorial Optimization, EGRES). |
|
8 * |
|
9 * Permission to use, modify and distribute this software is granted |
|
10 * provided that this copyright notice appears in all copies. For |
|
11 * precise terms see the accompanying LICENSE file. |
|
12 * |
|
13 * This software is provided "AS IS" with no warranty of any kind, |
|
14 * express or implied, and with no claim as to its suitability for any |
|
15 * purpose. |
|
16 * |
|
17 */ |
|
18 |
|
19 #ifndef LEMON_MIN_MEAN_CYCLE_H |
|
20 #define LEMON_MIN_MEAN_CYCLE_H |
|
21 |
|
22 /// \ingroup shortest_path |
|
23 /// |
|
24 /// \file |
|
25 /// \brief Howard's algorithm for finding a minimum mean cycle. |
|
26 |
|
27 #include <vector> |
|
28 #include <lemon/core.h> |
|
29 #include <lemon/path.h> |
|
30 #include <lemon/tolerance.h> |
|
31 #include <lemon/connectivity.h> |
|
32 |
|
33 namespace lemon { |
|
34 |
|
35 /// \addtogroup shortest_path |
|
36 /// @{ |
|
37 |
|
38 /// \brief Implementation of Howard's algorithm for finding a minimum |
|
39 /// mean cycle. |
|
40 /// |
|
41 /// \ref MinMeanCycle implements Howard's algorithm for finding a |
|
42 /// directed cycle of minimum mean length (cost) in a digraph. |
|
43 /// |
|
44 /// \tparam GR The type of the digraph the algorithm runs on. |
|
45 /// \tparam LEN The type of the length map. The default |
|
46 /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>". |
|
47 /// |
|
48 /// \warning \c LEN::Value must be convertible to \c double. |
|
49 #ifdef DOXYGEN |
|
50 template <typename GR, typename LEN> |
|
51 #else |
|
52 template < typename GR, |
|
53 typename LEN = typename GR::template ArcMap<int> > |
|
54 #endif |
|
55 class MinMeanCycle |
|
56 { |
|
57 public: |
|
58 |
|
59 /// The type of the digraph the algorithm runs on |
|
60 typedef GR Digraph; |
|
61 /// The type of the length map |
|
62 typedef LEN LengthMap; |
|
63 /// The type of the arc lengths |
|
64 typedef typename LengthMap::Value Value; |
|
65 /// The type of the paths |
|
66 typedef lemon::Path<Digraph> Path; |
|
67 |
|
68 private: |
|
69 |
|
70 TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); |
|
71 |
|
72 // The digraph the algorithm runs on |
|
73 const Digraph &_gr; |
|
74 // The length of the arcs |
|
75 const LengthMap &_length; |
|
76 |
|
77 // The total length of the found cycle |
|
78 Value _cycle_length; |
|
79 // The number of arcs on the found cycle |
|
80 int _cycle_size; |
|
81 // The found cycle |
|
82 Path *_cycle_path; |
|
83 |
|
84 bool _local_path; |
|
85 bool _cycle_found; |
|
86 Node _cycle_node; |
|
87 |
|
88 typename Digraph::template NodeMap<bool> _reached; |
|
89 typename Digraph::template NodeMap<double> _dist; |
|
90 typename Digraph::template NodeMap<Arc> _policy; |
|
91 |
|
92 typename Digraph::template NodeMap<int> _comp; |
|
93 int _comp_num; |
|
94 |
|
95 std::vector<Node> _nodes; |
|
96 std::vector<Arc> _arcs; |
|
97 Tolerance<double> _tol; |
|
98 |
|
99 public: |
|
100 |
|
101 /// \brief Constructor. |
|
102 /// |
|
103 /// The constructor of the class. |
|
104 /// |
|
105 /// \param digraph The digraph the algorithm runs on. |
|
106 /// \param length The lengths (costs) of the arcs. |
|
107 MinMeanCycle( const Digraph &digraph, |
|
108 const LengthMap &length ) : |
|
109 _gr(digraph), _length(length), _cycle_length(0), _cycle_size(-1), |
|
110 _cycle_path(NULL), _local_path(false), _reached(digraph), |
|
111 _dist(digraph), _policy(digraph), _comp(digraph) |
|
112 {} |
|
113 |
|
114 /// Destructor. |
|
115 ~MinMeanCycle() { |
|
116 if (_local_path) delete _cycle_path; |
|
117 } |
|
118 |
|
119 /// \brief Set the path structure for storing the found cycle. |
|
120 /// |
|
121 /// This function sets an external path structure for storing the |
|
122 /// found cycle. |
|
123 /// |
|
124 /// If you don't call this function before calling \ref run() or |
|
125 /// \ref init(), it will allocate a local \ref Path "path" |
|
126 /// structure. The destuctor deallocates this automatically |
|
127 /// allocated object, of course. |
|
128 /// |
|
129 /// \note The algorithm calls only the \ref lemon::Path::addBack() |
|
130 /// "addBack()" function of the given path structure. |
|
131 /// |
|
132 /// \return <tt>(*this)</tt> |
|
133 /// |
|
134 /// \sa cycle() |
|
135 MinMeanCycle& cyclePath(Path &path) { |
|
136 if (_local_path) { |
|
137 delete _cycle_path; |
|
138 _local_path = false; |
|
139 } |
|
140 _cycle_path = &path; |
|
141 return *this; |
|
142 } |
|
143 |
|
144 /// \name Execution control |
|
145 /// The simplest way to execute the algorithm is to call the \ref run() |
|
146 /// function.\n |
|
147 /// If you only need the minimum mean length, you may call \ref init() |
|
148 /// and \ref findMinMean(). |
|
149 /// If you would like to run the algorithm again (e.g. the underlying |
|
150 /// digraph and/or the arc lengths has been modified), you may not |
|
151 /// create a new instance of the class, rather call \ref reset(), |
|
152 /// \ref findMinMean() and \ref findCycle() instead. |
|
153 |
|
154 /// @{ |
|
155 |
|
156 /// \brief Run the algorithm. |
|
157 /// |
|
158 /// This function runs the algorithm. |
|
159 /// |
|
160 /// \return \c true if a directed cycle exists in the digraph. |
|
161 /// |
|
162 /// \note Apart from the return value, <tt>mmc.run()</tt> is just a |
|
163 /// shortcut of the following code. |
|
164 /// \code |
|
165 /// mmc.init(); |
|
166 /// mmc.findMinMean(); |
|
167 /// mmc.findCycle(); |
|
168 /// \endcode |
|
169 bool run() { |
|
170 init(); |
|
171 return findMinMean() && findCycle(); |
|
172 } |
|
173 |
|
174 /// \brief Initialize the internal data structures. |
|
175 /// |
|
176 /// This function initializes the internal data structures. |
|
177 /// |
|
178 /// \sa reset() |
|
179 void init() { |
|
180 _tol.epsilon(1e-6); |
|
181 if (!_cycle_path) { |
|
182 _local_path = true; |
|
183 _cycle_path = new Path; |
|
184 } |
|
185 _cycle_found = false; |
|
186 _comp_num = stronglyConnectedComponents(_gr, _comp); |
|
187 } |
|
188 |
|
189 /// \brief Reset the internal data structures. |
|
190 /// |
|
191 /// This function resets the internal data structures so that |
|
192 /// findMinMean() and findCycle() can be called again (e.g. when the |
|
193 /// underlying digraph and/or the arc lengths has been modified). |
|
194 /// |
|
195 /// \sa init() |
|
196 void reset() { |
|
197 if (_cycle_path) _cycle_path->clear(); |
|
198 _cycle_found = false; |
|
199 _comp_num = stronglyConnectedComponents(_gr, _comp); |
|
200 } |
|
201 |
|
202 /// \brief Find the minimum cycle mean. |
|
203 /// |
|
204 /// This function computes all the required data and finds the |
|
205 /// minimum mean length of the directed cycles in the digraph. |
|
206 /// |
|
207 /// \return \c true if a directed cycle exists in the digraph. |
|
208 /// |
|
209 /// \pre \ref init() must be called before using this function. |
|
210 bool findMinMean() { |
|
211 // Find the minimum cycle mean in the components |
|
212 for (int comp = 0; comp < _comp_num; ++comp) { |
|
213 if (!initCurrentComponent(comp)) continue; |
|
214 while (true) { |
|
215 if (!findPolicyCycles()) break; |
|
216 contractPolicyGraph(comp); |
|
217 if (!computeNodeDistances()) break; |
|
218 } |
|
219 } |
|
220 return _cycle_found; |
|
221 } |
|
222 |
|
223 /// \brief Find a minimum mean directed cycle. |
|
224 /// |
|
225 /// This function finds a directed cycle of minimum mean length |
|
226 /// in the digraph using the data computed by findMinMean(). |
|
227 /// |
|
228 /// \return \c true if a directed cycle exists in the digraph. |
|
229 /// |
|
230 /// \pre \ref init() and \ref findMinMean() must be called before |
|
231 /// using this function. |
|
232 bool findCycle() { |
|
233 if (!_cycle_found) return false; |
|
234 _cycle_path->addBack(_policy[_cycle_node]); |
|
235 for ( Node v = _cycle_node; |
|
236 (v = _gr.target(_policy[v])) != _cycle_node; ) { |
|
237 _cycle_path->addBack(_policy[v]); |
|
238 } |
|
239 return true; |
|
240 } |
|
241 |
|
242 /// @} |
|
243 |
|
244 /// \name Query Functions |
|
245 /// The result of the algorithm can be obtained using these |
|
246 /// functions.\n |
|
247 /// The algorithm should be executed before using them. |
|
248 |
|
249 /// @{ |
|
250 |
|
251 /// \brief Return the total length of the found cycle. |
|
252 /// |
|
253 /// This function returns the total length of the found cycle. |
|
254 /// |
|
255 /// \pre \ref run() or \ref findCycle() must be called before |
|
256 /// using this function. |
|
257 Value cycleLength() const { |
|
258 return _cycle_length; |
|
259 } |
|
260 |
|
261 /// \brief Return the number of arcs on the found cycle. |
|
262 /// |
|
263 /// This function returns the number of arcs on the found cycle. |
|
264 /// |
|
265 /// \pre \ref run() or \ref findCycle() must be called before |
|
266 /// using this function. |
|
267 int cycleArcNum() const { |
|
268 return _cycle_size; |
|
269 } |
|
270 |
|
271 /// \brief Return the mean length of the found cycle. |
|
272 /// |
|
273 /// This function returns the mean length of the found cycle. |
|
274 /// |
|
275 /// \note <tt>mmc.cycleMean()</tt> is just a shortcut of the |
|
276 /// following code. |
|
277 /// \code |
|
278 /// return double(mmc.cycleLength()) / mmc.cycleArcNum(); |
|
279 /// \endcode |
|
280 /// |
|
281 /// \pre \ref run() or \ref findMinMean() must be called before |
|
282 /// using this function. |
|
283 double cycleMean() const { |
|
284 return double(_cycle_length) / _cycle_size; |
|
285 } |
|
286 |
|
287 /// \brief Return the found cycle. |
|
288 /// |
|
289 /// This function returns a const reference to the path structure |
|
290 /// storing the found cycle. |
|
291 /// |
|
292 /// \pre \ref run() or \ref findCycle() must be called before using |
|
293 /// this function. |
|
294 /// |
|
295 /// \sa cyclePath() |
|
296 const Path& cycle() const { |
|
297 return *_cycle_path; |
|
298 } |
|
299 |
|
300 ///@} |
|
301 |
|
302 private: |
|
303 |
|
304 // Initialize the internal data structures for the current strongly |
|
305 // connected component and create the policy graph. |
|
306 // The policy graph can be represented by the _policy map because |
|
307 // the out-degree of every node is 1. |
|
308 bool initCurrentComponent(int comp) { |
|
309 // Find the nodes of the current component |
|
310 _nodes.clear(); |
|
311 for (NodeIt n(_gr); n != INVALID; ++n) { |
|
312 if (_comp[n] == comp) _nodes.push_back(n); |
|
313 } |
|
314 if (_nodes.size() <= 1) return false; |
|
315 // Find the arcs of the current component |
|
316 _arcs.clear(); |
|
317 for (ArcIt e(_gr); e != INVALID; ++e) { |
|
318 if ( _comp[_gr.source(e)] == comp && |
|
319 _comp[_gr.target(e)] == comp ) |
|
320 _arcs.push_back(e); |
|
321 } |
|
322 // Initialize _reached, _dist, _policy maps |
|
323 for (int i = 0; i < int(_nodes.size()); ++i) { |
|
324 _reached[_nodes[i]] = false; |
|
325 _policy[_nodes[i]] = INVALID; |
|
326 } |
|
327 Node u; Arc e; |
|
328 for (int j = 0; j < int(_arcs.size()); ++j) { |
|
329 e = _arcs[j]; |
|
330 u = _gr.source(e); |
|
331 if (!_reached[u] || _length[e] < _dist[u]) { |
|
332 _dist[u] = _length[e]; |
|
333 _policy[u] = e; |
|
334 _reached[u] = true; |
|
335 } |
|
336 } |
|
337 return true; |
|
338 } |
|
339 |
|
340 // Find all cycles in the policy graph. |
|
341 // Set _cycle_found to true if a cycle is found and set |
|
342 // _cycle_length, _cycle_size, _cycle_node to represent the minimum |
|
343 // mean cycle in the policy graph. |
|
344 bool findPolicyCycles() { |
|
345 typename Digraph::template NodeMap<int> level(_gr, -1); |
|
346 bool curr_cycle_found = false; |
|
347 Value clength; |
|
348 int csize; |
|
349 int path_cnt = 0; |
|
350 Node u, v; |
|
351 // Searching for cycles |
|
352 for (int i = 0; i < int(_nodes.size()); ++i) { |
|
353 if (level[_nodes[i]] < 0) { |
|
354 u = _nodes[i]; |
|
355 level[u] = path_cnt; |
|
356 while (level[u = _gr.target(_policy[u])] < 0) |
|
357 level[u] = path_cnt; |
|
358 if (level[u] == path_cnt) { |
|
359 // A cycle is found |
|
360 curr_cycle_found = true; |
|
361 clength = _length[_policy[u]]; |
|
362 csize = 1; |
|
363 for (v = u; (v = _gr.target(_policy[v])) != u; ) { |
|
364 clength += _length[_policy[v]]; |
|
365 ++csize; |
|
366 } |
|
367 if ( !_cycle_found || |
|
368 clength * _cycle_size < _cycle_length * csize ) { |
|
369 _cycle_found = true; |
|
370 _cycle_length = clength; |
|
371 _cycle_size = csize; |
|
372 _cycle_node = u; |
|
373 } |
|
374 } |
|
375 ++path_cnt; |
|
376 } |
|
377 } |
|
378 return curr_cycle_found; |
|
379 } |
|
380 |
|
381 // Contract the policy graph to be connected by cutting all cycles |
|
382 // except for the main cycle (i.e. the minimum mean cycle). |
|
383 void contractPolicyGraph(int comp) { |
|
384 // Find the component of the main cycle using reverse BFS search |
|
385 typename Digraph::template NodeMap<int> found(_gr, false); |
|
386 std::deque<Node> queue; |
|
387 queue.push_back(_cycle_node); |
|
388 found[_cycle_node] = true; |
|
389 Node u, v; |
|
390 while (!queue.empty()) { |
|
391 v = queue.front(); queue.pop_front(); |
|
392 for (InArcIt e(_gr, v); e != INVALID; ++e) { |
|
393 u = _gr.source(e); |
|
394 if (_policy[u] == e && !found[u]) { |
|
395 found[u] = true; |
|
396 queue.push_back(u); |
|
397 } |
|
398 } |
|
399 } |
|
400 // Connect all other nodes to this component using reverse BFS search |
|
401 queue.clear(); |
|
402 for (int i = 0; i < int(_nodes.size()); ++i) |
|
403 if (found[_nodes[i]]) queue.push_back(_nodes[i]); |
|
404 int found_cnt = queue.size(); |
|
405 while (found_cnt < int(_nodes.size())) { |
|
406 v = queue.front(); queue.pop_front(); |
|
407 for (InArcIt e(_gr, v); e != INVALID; ++e) { |
|
408 u = _gr.source(e); |
|
409 if (_comp[u] == comp && !found[u]) { |
|
410 found[u] = true; |
|
411 ++found_cnt; |
|
412 _policy[u] = e; |
|
413 queue.push_back(u); |
|
414 } |
|
415 } |
|
416 } |
|
417 } |
|
418 |
|
419 // Compute node distances in the policy graph and update the |
|
420 // policy graph if the node distances can be improved. |
|
421 bool computeNodeDistances() { |
|
422 // Compute node distances using reverse BFS search |
|
423 double cycle_mean = double(_cycle_length) / _cycle_size; |
|
424 typename Digraph::template NodeMap<int> found(_gr, false); |
|
425 std::deque<Node> queue; |
|
426 queue.push_back(_cycle_node); |
|
427 found[_cycle_node] = true; |
|
428 _dist[_cycle_node] = 0; |
|
429 Node u, v; |
|
430 while (!queue.empty()) { |
|
431 v = queue.front(); queue.pop_front(); |
|
432 for (InArcIt e(_gr, v); e != INVALID; ++e) { |
|
433 u = _gr.source(e); |
|
434 if (_policy[u] == e && !found[u]) { |
|
435 found[u] = true; |
|
436 _dist[u] = _dist[v] + _length[e] - cycle_mean; |
|
437 queue.push_back(u); |
|
438 } |
|
439 } |
|
440 } |
|
441 // Improving node distances |
|
442 bool improved = false; |
|
443 for (int j = 0; j < int(_arcs.size()); ++j) { |
|
444 Arc e = _arcs[j]; |
|
445 u = _gr.source(e); v = _gr.target(e); |
|
446 double delta = _dist[v] + _length[e] - cycle_mean; |
|
447 if (_tol.less(delta, _dist[u])) { |
|
448 improved = true; |
|
449 _dist[u] = delta; |
|
450 _policy[u] = e; |
|
451 } |
|
452 } |
|
453 return improved; |
|
454 } |
|
455 |
|
456 }; //class MinMeanCycle |
|
457 |
|
458 ///@} |
|
459 |
|
460 } //namespace lemon |
|
461 |
|
462 #endif //LEMON_MIN_MEAN_CYCLE_H |