148 TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); |
148 TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); |
149 |
149 |
150 // Data sturcture for path data |
150 // Data sturcture for path data |
151 struct PathData |
151 struct PathData |
152 { |
152 { |
153 bool found; |
|
154 LargeValue dist; |
153 LargeValue dist; |
155 Arc pred; |
154 Arc pred; |
156 PathData(bool f = false, LargeValue d = 0, Arc p = INVALID) : |
155 PathData(LargeValue d, Arc p = INVALID) : |
157 found(f), dist(d), pred(p) {} |
156 dist(d), pred(p) {} |
158 }; |
157 }; |
159 |
158 |
160 typedef typename Digraph::template NodeMap<std::vector<PathData> > |
159 typedef typename Digraph::template NodeMap<std::vector<PathData> > |
161 PathDataNodeMap; |
160 PathDataNodeMap; |
162 |
161 |
188 PathDataNodeMap _data; |
187 PathDataNodeMap _data; |
189 // The processed nodes in the last round |
188 // The processed nodes in the last round |
190 std::vector<Node> _process; |
189 std::vector<Node> _process; |
191 |
190 |
192 Tolerance _tolerance; |
191 Tolerance _tolerance; |
|
192 |
|
193 // Infinite constant |
|
194 const LargeValue INF; |
193 |
195 |
194 public: |
196 public: |
195 |
197 |
196 /// \name Named Template Parameters |
198 /// \name Named Template Parameters |
197 /// @{ |
199 /// @{ |
243 /// \param length The lengths (costs) of the arcs. |
245 /// \param length The lengths (costs) of the arcs. |
244 HartmannOrlin( const Digraph &digraph, |
246 HartmannOrlin( const Digraph &digraph, |
245 const LengthMap &length ) : |
247 const LengthMap &length ) : |
246 _gr(digraph), _length(length), _comp(digraph), _out_arcs(digraph), |
248 _gr(digraph), _length(length), _comp(digraph), _out_arcs(digraph), |
247 _best_found(false), _best_length(0), _best_size(1), |
249 _best_found(false), _best_length(0), _best_size(1), |
248 _cycle_path(NULL), _local_path(false), _data(digraph) |
250 _cycle_path(NULL), _local_path(false), _data(digraph), |
|
251 INF(std::numeric_limits<LargeValue>::has_infinity ? |
|
252 std::numeric_limits<LargeValue>::infinity() : |
|
253 std::numeric_limits<LargeValue>::max()) |
249 {} |
254 {} |
250 |
255 |
251 /// Destructor. |
256 /// Destructor. |
252 ~HartmannOrlin() { |
257 ~HartmannOrlin() { |
253 if (_local_path) delete _cycle_path; |
258 if (_local_path) delete _cycle_path; |
470 int n = _nodes->size(); |
475 int n = _nodes->size(); |
471 if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) { |
476 if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) { |
472 return false; |
477 return false; |
473 } |
478 } |
474 for (int i = 0; i < n; ++i) { |
479 for (int i = 0; i < n; ++i) { |
475 _data[(*_nodes)[i]].resize(n + 1); |
480 _data[(*_nodes)[i]].resize(n + 1, PathData(INF)); |
476 } |
481 } |
477 return true; |
482 return true; |
478 } |
483 } |
479 |
484 |
480 // Process all rounds of computing path data for the current component. |
485 // Process all rounds of computing path data for the current component. |
481 // _data[v][k] is the length of a shortest directed walk from the root |
486 // _data[v][k] is the length of a shortest directed walk from the root |
482 // node to node v containing exactly k arcs. |
487 // node to node v containing exactly k arcs. |
483 void processRounds() { |
488 void processRounds() { |
484 Node start = (*_nodes)[0]; |
489 Node start = (*_nodes)[0]; |
485 _data[start][0] = PathData(true, 0); |
490 _data[start][0] = PathData(0); |
486 _process.clear(); |
491 _process.clear(); |
487 _process.push_back(start); |
492 _process.push_back(start); |
488 |
493 |
489 int k, n = _nodes->size(); |
494 int k, n = _nodes->size(); |
490 int next_check = 4; |
495 int next_check = 4; |
515 u = _process[i]; |
520 u = _process[i]; |
516 for (int j = 0; j < int(_out_arcs[u].size()); ++j) { |
521 for (int j = 0; j < int(_out_arcs[u].size()); ++j) { |
517 e = _out_arcs[u][j]; |
522 e = _out_arcs[u][j]; |
518 v = _gr.target(e); |
523 v = _gr.target(e); |
519 d = _data[u][k-1].dist + _length[e]; |
524 d = _data[u][k-1].dist + _length[e]; |
520 if (!_data[v][k].found) { |
525 if (_tolerance.less(d, _data[v][k].dist)) { |
521 next.push_back(v); |
526 if (_data[v][k].dist == INF) next.push_back(v); |
522 _data[v][k] = PathData(true, _data[u][k-1].dist + _length[e], e); |
527 _data[v][k] = PathData(d, e); |
523 } |
|
524 else if (_tolerance.less(d, _data[v][k].dist)) { |
|
525 _data[v][k] = PathData(true, d, e); |
|
526 } |
528 } |
527 } |
529 } |
528 } |
530 } |
529 _process.swap(next); |
531 _process.swap(next); |
530 } |
532 } |
538 u = (*_nodes)[i]; |
540 u = (*_nodes)[i]; |
539 for (int j = 0; j < int(_out_arcs[u].size()); ++j) { |
541 for (int j = 0; j < int(_out_arcs[u].size()); ++j) { |
540 e = _out_arcs[u][j]; |
542 e = _out_arcs[u][j]; |
541 v = _gr.target(e); |
543 v = _gr.target(e); |
542 d = _data[u][k-1].dist + _length[e]; |
544 d = _data[u][k-1].dist + _length[e]; |
543 if (!_data[v][k].found || _tolerance.less(d, _data[v][k].dist)) { |
545 if (_tolerance.less(d, _data[v][k].dist)) { |
544 _data[v][k] = PathData(true, d, e); |
546 _data[v][k] = PathData(d, e); |
545 } |
547 } |
546 } |
548 } |
547 } |
549 } |
548 } |
550 } |
549 |
551 |
559 |
561 |
560 // Search for cycles that are already found |
562 // Search for cycles that are already found |
561 _curr_found = false; |
563 _curr_found = false; |
562 for (int i = 0; i < n; ++i) { |
564 for (int i = 0; i < n; ++i) { |
563 u = (*_nodes)[i]; |
565 u = (*_nodes)[i]; |
564 if (!_data[u][k].found) continue; |
566 if (_data[u][k].dist == INF) continue; |
565 for (int j = k; j >= 0; --j) { |
567 for (int j = k; j >= 0; --j) { |
566 if (level[u].first == i && level[u].second > 0) { |
568 if (level[u].first == i && level[u].second > 0) { |
567 // A cycle is found |
569 // A cycle is found |
568 length = _data[u][level[u].second].dist - _data[u][j].dist; |
570 length = _data[u][level[u].second].dist - _data[u][j].dist; |
569 size = level[u].second - j; |
571 size = level[u].second - j; |
584 LargeValue d; |
586 LargeValue d; |
585 if (_curr_found && k < n) { |
587 if (_curr_found && k < n) { |
586 // Find node potentials |
588 // Find node potentials |
587 for (int i = 0; i < n; ++i) { |
589 for (int i = 0; i < n; ++i) { |
588 u = (*_nodes)[i]; |
590 u = (*_nodes)[i]; |
589 pi[u] = std::numeric_limits<LargeValue>::max(); |
591 pi[u] = INF; |
590 for (int j = 0; j <= k; ++j) { |
592 for (int j = 0; j <= k; ++j) { |
591 d = _data[u][j].dist * _curr_size - j * _curr_length; |
593 if (_data[u][j].dist < INF) { |
592 if (_data[u][j].found && _tolerance.less(d, pi[u])) { |
594 d = _data[u][j].dist * _curr_size - j * _curr_length; |
593 pi[u] = d; |
595 if (_tolerance.less(d, pi[u])) pi[u] = d; |
594 } |
596 } |
595 } |
597 } |
596 } |
598 } |
597 |
599 |
598 // Check the optimality condition for all arcs |
600 // Check the optimality condition for all arcs |