1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
3 * This file is a part of LEMON, a generic C++ optimization library.
5 * Copyright (C) 2003-2009
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
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.
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
19 #ifndef LEMON_SUURBALLE_H
20 #define LEMON_SUURBALLE_H
22 ///\ingroup shortest_path
24 ///\brief An algorithm for finding arc-disjoint paths between two
25 /// nodes having minimum total length.
29 #include <lemon/bin_heap.h>
30 #include <lemon/path.h>
31 #include <lemon/list_graph.h>
32 #include <lemon/maps.h>
36 /// \addtogroup shortest_path
39 /// \brief Algorithm for finding arc-disjoint paths between two nodes
40 /// having minimum total length.
42 /// \ref lemon::Suurballe "Suurballe" implements an algorithm for
43 /// finding arc-disjoint paths having minimum total length (cost)
44 /// from a given source node to a given target node in a digraph.
46 /// Note that this problem is a special case of the \ref min_cost_flow
47 /// "minimum cost flow problem". This implementation is actually an
48 /// efficient specialized version of the \ref CapacityScaling
49 /// "successive shortest path" algorithm directly for this problem.
50 /// Therefore this class provides query functions for flow values and
51 /// node potentials (the dual solution) just like the minimum cost flow
54 /// \tparam GR The digraph type the algorithm runs on.
55 /// \tparam LEN The type of the length map.
56 /// The default value is <tt>GR::ArcMap<int></tt>.
58 /// \warning Length values should be \e non-negative.
60 /// \note For finding \e node-disjoint paths, this algorithm can be used
61 /// along with the \ref SplitNodes adaptor.
63 template <typename GR, typename LEN>
65 template < typename GR,
66 typename LEN = typename GR::template ArcMap<int> >
70 TEMPLATE_DIGRAPH_TYPEDEFS(GR);
72 typedef ConstMap<Arc, int> ConstArcMap;
73 typedef typename GR::template NodeMap<Arc> PredMap;
77 /// The type of the digraph the algorithm runs on.
79 /// The type of the length map.
80 typedef LEN LengthMap;
81 /// The type of the lengths.
82 typedef typename LengthMap::Value Length;
84 /// The type of the flow map.
85 typedef GR::ArcMap<int> FlowMap;
86 /// The type of the potential map.
87 typedef GR::NodeMap<Length> PotentialMap;
89 /// The type of the flow map.
90 typedef typename Digraph::template ArcMap<int> FlowMap;
91 /// The type of the potential map.
92 typedef typename Digraph::template NodeMap<Length> PotentialMap;
95 /// The type of the path structures.
96 typedef SimplePath<GR> Path;
100 // ResidualDijkstra is a special implementation of the
101 // Dijkstra algorithm for finding shortest paths in the
102 // residual network with respect to the reduced arc lengths
103 // and modifying the node potentials according to the
104 // distance of the nodes.
105 class ResidualDijkstra
107 typedef typename Digraph::template NodeMap<int> HeapCrossRef;
108 typedef BinHeap<Length, HeapCrossRef> Heap;
112 const Digraph &_graph;
113 const LengthMap &_length;
114 const FlowMap &_flow;
121 std::vector<Node> _proc_nodes;
126 ResidualDijkstra(Suurballe &srb) :
127 _graph(srb._graph), _length(srb._length),
128 _flow(*srb._flow), _pi(*srb._potential), _pred(srb._pred),
129 _s(srb._s), _t(srb._t), _dist(_graph) {}
131 // Run the algorithm and return true if a path is found
132 // from the source node to the target node.
134 return cnt == 0 ? startFirst() : start();
139 // Execute the algorithm for the first time (the flow and potential
140 // functions have to be identically zero).
142 HeapCrossRef heap_cross_ref(_graph, Heap::PRE_HEAP);
143 Heap heap(heap_cross_ref);
149 while (!heap.empty() && heap.top() != _t) {
150 Node u = heap.top(), v;
151 Length d = heap.prio(), dn;
152 _dist[u] = heap.prio();
153 _proc_nodes.push_back(u);
156 // Traverse outgoing arcs
157 for (OutArcIt e(_graph, u); e != INVALID; ++e) {
158 v = _graph.target(e);
159 switch(heap.state(v)) {
161 heap.push(v, d + _length[e]);
167 heap.decrease(v, dn);
171 case Heap::POST_HEAP:
176 if (heap.empty()) return false;
178 // Update potentials of processed nodes
179 Length t_dist = heap.prio();
180 for (int i = 0; i < int(_proc_nodes.size()); ++i)
181 _pi[_proc_nodes[i]] = _dist[_proc_nodes[i]] - t_dist;
185 // Execute the algorithm.
187 HeapCrossRef heap_cross_ref(_graph, Heap::PRE_HEAP);
188 Heap heap(heap_cross_ref);
194 while (!heap.empty() && heap.top() != _t) {
195 Node u = heap.top(), v;
196 Length d = heap.prio() + _pi[u], dn;
197 _dist[u] = heap.prio();
198 _proc_nodes.push_back(u);
201 // Traverse outgoing arcs
202 for (OutArcIt e(_graph, u); e != INVALID; ++e) {
204 v = _graph.target(e);
205 switch(heap.state(v)) {
207 heap.push(v, d + _length[e] - _pi[v]);
211 dn = d + _length[e] - _pi[v];
213 heap.decrease(v, dn);
217 case Heap::POST_HEAP:
223 // Traverse incoming arcs
224 for (InArcIt e(_graph, u); e != INVALID; ++e) {
226 v = _graph.source(e);
227 switch(heap.state(v)) {
229 heap.push(v, d - _length[e] - _pi[v]);
233 dn = d - _length[e] - _pi[v];
235 heap.decrease(v, dn);
239 case Heap::POST_HEAP:
245 if (heap.empty()) return false;
247 // Update potentials of processed nodes
248 Length t_dist = heap.prio();
249 for (int i = 0; i < int(_proc_nodes.size()); ++i)
250 _pi[_proc_nodes[i]] += _dist[_proc_nodes[i]] - t_dist;
254 }; //class ResidualDijkstra
258 // The digraph the algorithm runs on
259 const Digraph &_graph;
261 const LengthMap &_length;
263 // Arc map of the current flow
266 // Node map of the current potentials
267 PotentialMap *_potential;
268 bool _local_potential;
275 // Container to store the found paths
276 std::vector<Path> _paths;
284 /// \brief Constructor.
288 /// \param graph The digraph the algorithm runs on.
289 /// \param length The length (cost) values of the arcs.
290 Suurballe( const Digraph &graph,
291 const LengthMap &length ) :
292 _graph(graph), _length(length), _flow(0), _local_flow(false),
293 _potential(0), _local_potential(false), _pred(graph)
298 if (_local_flow) delete _flow;
299 if (_local_potential) delete _potential;
302 /// \brief Set the flow map.
304 /// This function sets the flow map.
305 /// If it is not used before calling \ref run() or \ref init(),
306 /// an instance will be allocated automatically. The destructor
307 /// deallocates this automatically allocated map, of course.
309 /// The found flow contains only 0 and 1 values, since it is the
310 /// union of the found arc-disjoint paths.
312 /// \return <tt>(*this)</tt>
313 Suurballe& flowMap(FlowMap &map) {
322 /// \brief Set the potential map.
324 /// This function sets the potential map.
325 /// If it is not used before calling \ref run() or \ref init(),
326 /// an instance will be allocated automatically. The destructor
327 /// deallocates this automatically allocated map, of course.
329 /// The node potentials provide the dual solution of the underlying
330 /// \ref min_cost_flow "minimum cost flow problem".
332 /// \return <tt>(*this)</tt>
333 Suurballe& potentialMap(PotentialMap &map) {
334 if (_local_potential) {
336 _local_potential = false;
342 /// \name Execution Control
343 /// The simplest way to execute the algorithm is to call the run()
346 /// If you only need the flow that is the union of the found
347 /// arc-disjoint paths, you may call init() and findFlow().
351 /// \brief Run the algorithm.
353 /// This function runs the algorithm.
355 /// \param s The source node.
356 /// \param t The target node.
357 /// \param k The number of paths to be found.
359 /// \return \c k if there are at least \c k arc-disjoint paths from
360 /// \c s to \c t in the digraph. Otherwise it returns the number of
361 /// arc-disjoint paths found.
363 /// \note Apart from the return value, <tt>s.run(s, t, k)</tt> is
364 /// just a shortcut of the following code.
367 /// s.findFlow(t, k);
370 int run(const Node& s, const Node& t, int k = 2) {
377 /// \brief Initialize the algorithm.
379 /// This function initializes the algorithm.
381 /// \param s The source node.
382 void init(const Node& s) {
387 _flow = new FlowMap(_graph);
391 _potential = new PotentialMap(_graph);
392 _local_potential = true;
394 for (ArcIt e(_graph); e != INVALID; ++e) (*_flow)[e] = 0;
395 for (NodeIt n(_graph); n != INVALID; ++n) (*_potential)[n] = 0;
398 /// \brief Execute the algorithm to find an optimal flow.
400 /// This function executes the successive shortest path algorithm to
401 /// find a minimum cost flow, which is the union of \c k (or less)
402 /// arc-disjoint paths.
404 /// \param t The target node.
405 /// \param k The number of paths to be found.
407 /// \return \c k if there are at least \c k arc-disjoint paths from
408 /// the source node to the given node \c t in the digraph.
409 /// Otherwise it returns the number of arc-disjoint paths found.
411 /// \pre \ref init() must be called before using this function.
412 int findFlow(const Node& t, int k = 2) {
414 ResidualDijkstra dijkstra(*this);
416 // Find shortest paths
418 while (_path_num < k) {
420 if (!dijkstra.run(_path_num)) break;
423 // Set the flow along the found shortest path
426 while ((e = _pred[u]) != INVALID) {
427 if (u == _graph.target(e)) {
429 u = _graph.source(e);
432 u = _graph.target(e);
439 /// \brief Compute the paths from the flow.
441 /// This function computes arc-disjoint paths from the found minimum
442 /// cost flow, which is the union of them.
444 /// \pre \ref init() and \ref findFlow() must be called before using
447 FlowMap res_flow(_graph);
448 for(ArcIt a(_graph); a != INVALID; ++a) res_flow[a] = (*_flow)[a];
451 _paths.resize(_path_num);
452 for (int i = 0; i < _path_num; ++i) {
455 OutArcIt e(_graph, n);
456 for ( ; res_flow[e] == 0; ++e) ;
457 n = _graph.target(e);
458 _paths[i].addBack(e);
466 /// \name Query Functions
467 /// The results of the algorithm can be obtained using these
469 /// \n The algorithm should be executed before using them.
473 /// \brief Return the total length of the found paths.
475 /// This function returns the total length of the found paths, i.e.
476 /// the total cost of the found flow.
477 /// The complexity of the function is O(e).
479 /// \pre \ref run() or \ref findFlow() must be called before using
481 Length totalLength() const {
483 for (ArcIt e(_graph); e != INVALID; ++e)
484 c += (*_flow)[e] * _length[e];
488 /// \brief Return the flow value on the given arc.
490 /// This function returns the flow value on the given arc.
491 /// It is \c 1 if the arc is involved in one of the found arc-disjoint
492 /// paths, otherwise it is \c 0.
494 /// \pre \ref run() or \ref findFlow() must be called before using
496 int flow(const Arc& arc) const {
497 return (*_flow)[arc];
500 /// \brief Return a const reference to an arc map storing the
503 /// This function returns a const reference to an arc map storing
504 /// the flow that is the union of the found arc-disjoint paths.
506 /// \pre \ref run() or \ref findFlow() must be called before using
508 const FlowMap& flowMap() const {
512 /// \brief Return the potential of the given node.
514 /// This function returns the potential of the given node.
515 /// The node potentials provide the dual solution of the
516 /// underlying \ref min_cost_flow "minimum cost flow problem".
518 /// \pre \ref run() or \ref findFlow() must be called before using
520 Length potential(const Node& node) const {
521 return (*_potential)[node];
524 /// \brief Return a const reference to a node map storing the
525 /// found potentials (the dual solution).
527 /// This function returns a const reference to a node map storing
528 /// the found potentials that provide the dual solution of the
529 /// underlying \ref min_cost_flow "minimum cost flow problem".
531 /// \pre \ref run() or \ref findFlow() must be called before using
533 const PotentialMap& potentialMap() const {
537 /// \brief Return the number of the found paths.
539 /// This function returns the number of the found paths.
541 /// \pre \ref run() or \ref findFlow() must be called before using
543 int pathNum() const {
547 /// \brief Return a const reference to the specified path.
549 /// This function returns a const reference to the specified path.
551 /// \param i The function returns the <tt>i</tt>-th path.
552 /// \c i must be between \c 0 and <tt>%pathNum()-1</tt>.
554 /// \pre \ref run() or \ref findPaths() must be called before using
556 const Path& path(int i) const {
568 #endif //LEMON_SUURBALLE_H