Revert to long long int since currently I don't know a better solution.
3 * This file is a part of LEMON, a generic C++ optimization library
5 * Copyright (C) 2003-2008
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 edge-disjoint paths between two
25 /// nodes having minimum total length.
28 #include <lemon/bin_heap.h>
29 #include <lemon/path.h>
33 /// \addtogroup shortest_path
36 /// \brief Implementation of an algorithm for finding edge-disjoint
37 /// paths between two nodes having minimum total length.
39 /// \ref lemon::Suurballe "Suurballe" implements an algorithm for
40 /// finding edge-disjoint paths having minimum total length (cost)
41 /// from a given source node to a given target node in a directed
44 /// In fact, this implementation is the specialization of the
45 /// \ref CapacityScaling "successive shortest path" algorithm.
47 /// \tparam Graph The directed graph type the algorithm runs on.
48 /// \tparam LengthMap The type of the length (cost) map.
50 /// \warning Length values should be \e non-negative \e integers.
52 /// \note For finding node-disjoint paths this algorithm can be used
53 /// with \ref SplitGraphAdaptor.
55 /// \author Attila Bernath and Peter Kovacs
57 template < typename Graph,
58 typename LengthMap = typename Graph::template EdgeMap<int> >
61 GRAPH_TYPEDEFS(typename Graph);
63 typedef typename LengthMap::Value Length;
64 typedef ConstMap<Edge, int> ConstEdgeMap;
65 typedef typename Graph::template NodeMap<Edge> PredMap;
69 /// The type of the flow map.
70 typedef typename Graph::template EdgeMap<int> FlowMap;
71 /// The type of the potential map.
72 typedef typename Graph::template NodeMap<Length> PotentialMap;
73 /// The type of the path structures.
74 typedef SimplePath<Graph> Path;
78 /// \brief Special implementation of the \ref Dijkstra algorithm
79 /// for finding shortest paths in the residual network.
81 /// \ref ResidualDijkstra is a special implementation of the
82 /// \ref Dijkstra algorithm for finding shortest paths in the
83 /// residual network of the graph with respect to the reduced edge
84 /// lengths and modifying the node potentials according to the
85 /// distance of the nodes.
86 class ResidualDijkstra
88 typedef typename Graph::template NodeMap<int> HeapCrossRef;
89 typedef BinHeap<Length, HeapCrossRef> Heap;
93 // The directed graph the algorithm runs on
98 const LengthMap &_length;
99 PotentialMap &_potential;
105 // The processed (i.e. permanently labeled) nodes
106 std::vector<Node> _proc_nodes;
114 ResidualDijkstra( const Graph &graph,
116 const LengthMap &length,
117 PotentialMap &potential,
120 _graph(graph), _flow(flow), _length(length), _potential(potential),
121 _dist(graph), _pred(pred), _s(s), _t(t) {}
123 /// \brief Runs the algorithm. Returns \c true if a path is found
124 /// from the source node to the target node.
126 HeapCrossRef heap_cross_ref(_graph, Heap::PRE_HEAP);
127 Heap heap(heap_cross_ref);
133 while (!heap.empty() && heap.top() != _t) {
134 Node u = heap.top(), v;
135 Length d = heap.prio() + _potential[u], nd;
136 _dist[u] = heap.prio();
138 _proc_nodes.push_back(u);
140 // Traversing outgoing edges
141 for (OutEdgeIt e(_graph, u); e != INVALID; ++e) {
143 v = _graph.target(e);
144 switch(heap.state(v)) {
146 heap.push(v, d + _length[e] - _potential[v]);
150 nd = d + _length[e] - _potential[v];
152 heap.decrease(v, nd);
156 case Heap::POST_HEAP:
162 // Traversing incoming edges
163 for (InEdgeIt e(_graph, u); e != INVALID; ++e) {
165 v = _graph.source(e);
166 switch(heap.state(v)) {
168 heap.push(v, d - _length[e] - _potential[v]);
172 nd = d - _length[e] - _potential[v];
174 heap.decrease(v, nd);
178 case Heap::POST_HEAP:
184 if (heap.empty()) return false;
186 // Updating potentials of processed nodes
187 Length t_dist = heap.prio();
188 for (int i = 0; i < int(_proc_nodes.size()); ++i)
189 _potential[_proc_nodes[i]] += _dist[_proc_nodes[i]] - t_dist;
193 }; //class ResidualDijkstra
197 // The directed graph the algorithm runs on
200 const LengthMap &_length;
202 // Edge map of the current flow
205 // Node map of the current potentials
206 PotentialMap *_potential;
207 bool _local_potential;
214 // Container to store the found paths
215 std::vector< SimplePath<Graph> > paths;
220 // Implementation of the Dijkstra algorithm for finding augmenting
221 // shortest paths in the residual network
222 ResidualDijkstra *_dijkstra;
226 /// \brief Constructor.
230 /// \param graph The directed graph the algorithm runs on.
231 /// \param length The length (cost) values of the edges.
232 /// \param s The source node.
233 /// \param t The target node.
234 Suurballe( const Graph &graph,
235 const LengthMap &length,
237 _graph(graph), _length(length), _flow(0), _local_flow(false),
238 _potential(0), _local_potential(false), _source(s), _target(t),
243 if (_local_flow) delete _flow;
244 if (_local_potential) delete _potential;
248 /// \brief Sets the flow map.
250 /// Sets the flow map.
252 /// The found flow contains only 0 and 1 values. It is the union of
253 /// the found edge-disjoint paths.
255 /// \return \c (*this)
256 Suurballe& flowMap(FlowMap &map) {
265 /// \brief Sets the potential map.
267 /// Sets the potential map.
269 /// The potentials provide the dual solution of the underlying
270 /// minimum cost flow problem.
272 /// \return \c (*this)
273 Suurballe& potentialMap(PotentialMap &map) {
274 if (_local_potential) {
276 _local_potential = false;
282 /// \name Execution control
283 /// The simplest way to execute the algorithm is to call the run()
286 /// If you only need the flow that is the union of the found
287 /// edge-disjoint paths, you may call init() and findFlow().
291 /// \brief Runs the algorithm.
293 /// Runs the algorithm.
295 /// \param k The number of paths to be found.
297 /// \return \c k if there are at least \c k edge-disjoint paths
298 /// from \c s to \c t. Otherwise it returns the number of
299 /// edge-disjoint paths found.
301 /// \note Apart from the return value, <tt>s.run(k)</tt> is just a
302 /// shortcut of the following code.
315 /// \brief Initializes the algorithm.
317 /// Initializes the algorithm.
321 _flow = new FlowMap(_graph);
325 _potential = new PotentialMap(_graph);
326 _local_potential = true;
328 for (EdgeIt e(_graph); e != INVALID; ++e) (*_flow)[e] = 0;
329 for (NodeIt n(_graph); n != INVALID; ++n) (*_potential)[n] = 0;
331 _dijkstra = new ResidualDijkstra( _graph, *_flow, _length,
336 /// \brief Executes the successive shortest path algorithm to find
339 /// Executes the successive shortest path algorithm to find a
340 /// minimum cost flow, which is the union of \c k or less
341 /// edge-disjoint paths.
343 /// \return \c k if there are at least \c k edge-disjoint paths
344 /// from \c s to \c t. Otherwise it returns the number of
345 /// edge-disjoint paths found.
347 /// \pre \ref init() must be called before using this function.
348 int findFlow(int k = 2) {
349 // Finding shortest paths
351 while (_path_num < k) {
353 if (!_dijkstra->run()) break;
356 // Setting the flow along the found shortest path
359 while ((e = _pred[u]) != INVALID) {
360 if (u == _graph.target(e)) {
362 u = _graph.source(e);
365 u = _graph.target(e);
372 /// \brief Computes the paths from the flow.
374 /// Computes the paths from the flow.
376 /// \pre \ref init() and \ref findFlow() must be called before using
379 // Creating the residual flow map (the union of the paths not
381 FlowMap res_flow(*_flow);
384 paths.resize(_path_num);
385 for (int i = 0; i < _path_num; ++i) {
387 while (n != _target) {
388 OutEdgeIt e(_graph, n);
389 for ( ; res_flow[e] == 0; ++e) ;
390 n = _graph.target(e);
399 /// \name Query Functions
400 /// The result of the algorithm can be obtained using these
402 /// \n The algorithm should be executed before using them.
406 /// \brief Returns a const reference to the edge map storing the
409 /// Returns a const reference to the edge map storing the flow that
410 /// is the union of the found edge-disjoint paths.
412 /// \pre \ref run() or findFlow() must be called before using this
414 const FlowMap& flowMap() const {
418 /// \brief Returns a const reference to the node map storing the
419 /// found potentials (the dual solution).
421 /// Returns a const reference to the node map storing the found
422 /// potentials that provide the dual solution of the underlying
423 /// minimum cost flow problem.
425 /// \pre \ref run() or findFlow() must be called before using this
427 const PotentialMap& potentialMap() const {
431 /// \brief Returns the flow on the given edge.
433 /// Returns the flow on the given edge.
434 /// It is \c 1 if the edge is involved in one of the found paths,
435 /// otherwise it is \c 0.
437 /// \pre \ref run() or findFlow() must be called before using this
439 int flow(const Edge& edge) const {
440 return (*_flow)[edge];
443 /// \brief Returns the potential of the given node.
445 /// Returns the potential of the given node.
447 /// \pre \ref run() or findFlow() must be called before using this
449 Length potential(const Node& node) const {
450 return (*_potential)[node];
453 /// \brief Returns the total length (cost) of the found paths (flow).
455 /// Returns the total length (cost) of the found paths (flow).
456 /// The complexity of the function is \f$ O(e) \f$.
458 /// \pre \ref run() or findFlow() must be called before using this
460 Length totalLength() const {
462 for (EdgeIt e(_graph); e != INVALID; ++e)
463 c += (*_flow)[e] * _length[e];
467 /// \brief Returns the number of the found paths.
469 /// Returns the number of the found paths.
471 /// \pre \ref run() or findFlow() must be called before using this
473 int pathNum() const {
477 /// \brief Returns a const reference to the specified path.
479 /// Returns a const reference to the specified path.
481 /// \param i The function returns the \c i-th path.
482 /// \c i must be between \c 0 and <tt>%pathNum()-1</tt>.
484 /// \pre \ref run() or findPaths() must be called before using this
486 Path path(int i) const {
498 #endif //LEMON_SUURBALLE_H