COIN-OR::LEMON - Graph Library

source: lemon-0.x/lemon/edmonds_karp.h @ 2553:bfced05fa852

Last change on this file since 2553:bfced05fa852 was 2553:bfced05fa852, checked in by Alpar Juttner, 12 years ago

Happy New Year to LEMON (+ better update-copyright-header script)

File size: 14.6 KB
Line 
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_EDMONDS_KARP_H
20#define LEMON_EDMONDS_KARP_H
21
22/// \file
23/// \ingroup max_flow
24/// \brief Implementation of the Edmonds-Karp algorithm.
25
26#include <lemon/tolerance.h>
27#include <vector>
28
29namespace lemon {
30
31  /// \brief Default traits class of EdmondsKarp class.
32  ///
33  /// Default traits class of EdmondsKarp class.
34  /// \param _Graph Graph type.
35  /// \param _CapacityMap Type of capacity map.
36  template <typename _Graph, typename _CapacityMap>
37  struct EdmondsKarpDefaultTraits {
38
39    /// \brief The graph type the algorithm runs on.
40    typedef _Graph Graph;
41
42    /// \brief The type of the map that stores the edge capacities.
43    ///
44    /// The type of the map that stores the edge capacities.
45    /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
46    typedef _CapacityMap CapacityMap;
47
48    /// \brief The type of the length of the edges.
49    typedef typename CapacityMap::Value Value;
50
51    /// \brief The map type that stores the flow values.
52    ///
53    /// The map type that stores the flow values.
54    /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept.
55    typedef typename Graph::template EdgeMap<Value> FlowMap;
56
57    /// \brief Instantiates a FlowMap.
58    ///
59    /// This function instantiates a \ref FlowMap.
60    /// \param graph The graph, to which we would like to define the flow map.
61    static FlowMap* createFlowMap(const Graph& graph) {
62      return new FlowMap(graph);
63    }
64
65    /// \brief The tolerance used by the algorithm
66    ///
67    /// The tolerance used by the algorithm to handle inexact computation.
68    typedef Tolerance<Value> Tolerance;
69
70  };
71
72  /// \ingroup max_flow
73  ///
74  /// \brief Edmonds-Karp algorithms class.
75  ///
76  /// This class provides an implementation of the \e Edmonds-Karp \e
77  /// algorithm producing a flow of maximum value in a directed
78  /// graphs. The Edmonds-Karp algorithm is slower than the Preflow
79  /// algorithm but it has an advantage of the step-by-step execution
80  /// control with feasible flow solutions. The \e source node, the \e
81  /// target node, the \e capacity of the edges and the \e starting \e
82  /// flow value of the edges should be passed to the algorithm
83  /// through the constructor.
84  ///
85  /// The time complexity of the algorithm is \f$ O(nm^2) \f$ in
86  /// worst case.  Always try the preflow algorithm instead of this if
87  /// you just want to compute the optimal flow.
88  ///
89  /// \param _Graph The directed graph type the algorithm runs on.
90  /// \param _CapacityMap The capacity map type.
91  /// \param _Traits Traits class to set various data types used by
92  /// the algorithm.  The default traits class is \ref
93  /// EdmondsKarpDefaultTraits.  See \ref EdmondsKarpDefaultTraits for the
94  /// documentation of a Edmonds-Karp traits class.
95  ///
96  /// \author Balazs Dezso
97#ifdef DOXYGEN
98  template <typename _Graph, typename _CapacityMap, typename _Traits>
99#else
100  template <typename _Graph,
101            typename _CapacityMap = typename _Graph::template EdgeMap<int>,
102            typename _Traits = EdmondsKarpDefaultTraits<_Graph, _CapacityMap> >
103#endif
104  class EdmondsKarp {
105  public:
106
107    typedef _Traits Traits;
108    typedef typename Traits::Graph Graph;
109    typedef typename Traits::CapacityMap CapacityMap;
110    typedef typename Traits::Value Value;
111
112    typedef typename Traits::FlowMap FlowMap;
113    typedef typename Traits::Tolerance Tolerance;
114
115    /// \brief \ref Exception for the case when the source equals the target.
116    ///
117    /// \ref Exception for the case when the source equals the target.
118    ///
119    class InvalidArgument : public lemon::LogicError {
120    public:
121      virtual const char* what() const throw() {
122        return "lemon::EdmondsKarp::InvalidArgument";
123      }
124    };
125
126
127  private:
128
129    GRAPH_TYPEDEFS(typename Graph);
130    typedef typename Graph::template NodeMap<Edge> PredMap;
131   
132    const Graph& _graph;
133    const CapacityMap* _capacity;
134
135    Node _source, _target;
136
137    FlowMap* _flow;
138    bool _local_flow;
139
140    PredMap* _pred;
141    std::vector<Node> _queue;
142   
143    Tolerance _tolerance;
144    Value _flow_value;
145
146    void createStructures() {
147      if (!_flow) {
148        _flow = Traits::createFlowMap(_graph);
149        _local_flow = true;
150      }
151      if (!_pred) {
152        _pred = new PredMap(_graph);
153      }
154      _queue.resize(countNodes(_graph));
155    }
156
157    void destroyStructures() {
158      if (_local_flow) {
159        delete _flow;
160      }
161      if (_pred) {
162        delete _pred;
163      }
164    }
165   
166  public:
167
168    ///\name Named template parameters
169
170    ///@{
171
172    template <typename _FlowMap>
173    struct DefFlowMapTraits : public Traits {
174      typedef _FlowMap FlowMap;
175      static FlowMap *createFlowMap(const Graph&) {
176        throw UninitializedParameter();
177      }
178    };
179
180    /// \brief \ref named-templ-param "Named parameter" for setting
181    /// FlowMap type
182    ///
183    /// \ref named-templ-param "Named parameter" for setting FlowMap
184    /// type
185    template <typename _FlowMap>
186    struct DefFlowMap
187      : public EdmondsKarp<Graph, CapacityMap, DefFlowMapTraits<_FlowMap> > {
188      typedef EdmondsKarp<Graph, CapacityMap, DefFlowMapTraits<_FlowMap> >
189      Create;
190    };
191
192
193    /// @}
194
195  protected:
196   
197    EdmondsKarp() {}
198
199  public:
200
201    /// \brief The constructor of the class.
202    ///
203    /// The constructor of the class.
204    /// \param graph The directed graph the algorithm runs on.
205    /// \param capacity The capacity of the edges.
206    /// \param source The source node.
207    /// \param target The target node.
208    EdmondsKarp(const Graph& graph, const CapacityMap& capacity,
209                Node source, Node target)
210      : _graph(graph), _capacity(&capacity), _source(source), _target(target),
211        _flow(0), _local_flow(false), _pred(0), _tolerance(), _flow_value()
212    {
213      if (_source == _target) {
214        throw InvalidArgument();
215      }
216    }
217
218    /// \brief Destrcutor.
219    ///
220    /// Destructor.
221    ~EdmondsKarp() {
222      destroyStructures();
223    }
224
225    /// \brief Sets the capacity map.
226    ///
227    /// Sets the capacity map.
228    /// \return \c (*this)
229    EdmondsKarp& capacityMap(const CapacityMap& map) {
230      _capacity = &map;
231      return *this;
232    }
233
234    /// \brief Sets the flow map.
235    ///
236    /// Sets the flow map.
237    /// \return \c (*this)
238    EdmondsKarp& flowMap(FlowMap& map) {
239      if (_local_flow) {
240        delete _flow;
241        _local_flow = false;
242      }
243      _flow = &map;
244      return *this;
245    }
246
247    /// \brief Returns the flow map.
248    ///
249    /// \return The flow map.
250    const FlowMap& flowMap() {
251      return *_flow;
252    }
253
254    /// \brief Sets the source node.
255    ///
256    /// Sets the source node.
257    /// \return \c (*this)
258    EdmondsKarp& source(const Node& node) {
259      _source = node;
260      return *this;
261    }
262
263    /// \brief Sets the target node.
264    ///
265    /// Sets the target node.
266    /// \return \c (*this)
267    EdmondsKarp& target(const Node& node) {
268      _target = node;
269      return *this;
270    }
271
272    /// \brief Sets the tolerance used by algorithm.
273    ///
274    /// Sets the tolerance used by algorithm.
275    EdmondsKarp& tolerance(const Tolerance& tolerance) const {
276      _tolerance = tolerance;
277      return *this;
278    }
279
280    /// \brief Returns the tolerance used by algorithm.
281    ///
282    /// Returns the tolerance used by algorithm.
283    const Tolerance& tolerance() const {
284      return tolerance;
285    }
286
287    /// \name Execution control The simplest way to execute the
288    /// algorithm is to use the \c run() member functions.
289    /// \n
290    /// If you need more control on initial solution or
291    /// execution then you have to call one \ref init() function and then
292    /// the start() or multiple times the \c augment() member function. 
293   
294    ///@{
295
296    /// \brief Initializes the algorithm
297    ///
298    /// It sets the flow to empty flow.
299    void init() {
300      createStructures();
301      for (EdgeIt it(_graph); it != INVALID; ++it) {
302        _flow->set(it, 0);
303      }
304      _flow_value = 0;
305    }
306   
307    /// \brief Initializes the algorithm
308    ///
309    /// Initializes the flow to the \c flowMap. The \c flowMap should
310    /// contain a feasible flow, ie. in each node excluding the source
311    /// and the target the incoming flow should be equal to the
312    /// outgoing flow.
313    template <typename FlowMap>
314    void flowInit(const FlowMap& flowMap) {
315      createStructures();
316      for (EdgeIt e(_graph); e != INVALID; ++e) {
317        _flow->set(e, flowMap[e]);
318      }
319      _flow_value = 0;
320      for (OutEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
321        _flow_value += (*_flow)[jt];
322      }
323      for (InEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
324        _flow_value -= (*_flow)[jt];
325      }
326    }
327
328    /// \brief Initializes the algorithm
329    ///
330    /// Initializes the flow to the \c flowMap. The \c flowMap should
331    /// contain a feasible flow, ie. in each node excluding the source
332    /// and the target the incoming flow should be equal to the
333    /// outgoing flow. 
334    /// \return %False when the given flowMap does not contain
335    /// feasible flow.
336    template <typename FlowMap>
337    bool checkedFlowInit(const FlowMap& flowMap) {
338      createStructures();
339      for (EdgeIt e(_graph); e != INVALID; ++e) {
340        _flow->set(e, flowMap[e]);
341      }
342      for (NodeIt it(_graph); it != INVALID; ++it) {
343        if (it == _source || it == _target) continue;
344        Value outFlow = 0;
345        for (OutEdgeIt jt(_graph, it); jt != INVALID; ++jt) {
346          outFlow += (*_flow)[jt];
347        }
348        Value inFlow = 0;
349        for (InEdgeIt jt(_graph, it); jt != INVALID; ++jt) {
350          inFlow += (*_flow)[jt];
351        }
352        if (_tolerance.different(outFlow, inFlow)) {
353          return false;
354        }
355      }
356      for (EdgeIt it(_graph); it != INVALID; ++it) {
357        if (_tolerance.less((*_flow)[it], 0)) return false;
358        if (_tolerance.less((*_capacity)[it], (*_flow)[it])) return false;
359      }
360      _flow_value = 0;
361      for (OutEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
362        _flow_value += (*_flow)[jt];
363      }
364      for (InEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
365        _flow_value -= (*_flow)[jt];
366      }
367      return true;
368    }
369
370    /// \brief Augment the solution on an edge shortest path.
371    ///
372    /// Augment the solution on an edge shortest path. It search an
373    /// edge shortest path between the source and the target
374    /// in the residual graph with the bfs algoritm.
375    /// Then it increase the flow on this path with the minimal residual
376    /// capacity on the path. If there is not such path it gives back
377    /// false.
378    /// \return %False when the augmenting is not success so the
379    /// current flow is a feasible and optimal solution.
380    bool augment() {
381      for (NodeIt n(_graph); n != INVALID; ++n) {
382        _pred->set(n, INVALID);
383      }
384     
385      int first = 0, last = 1;
386     
387      _queue[0] = _source;
388      _pred->set(_source, OutEdgeIt(_graph, _source));
389
390      while (first != last && (*_pred)[_target] == INVALID) {
391        Node n = _queue[first++];
392       
393        for (OutEdgeIt e(_graph, n); e != INVALID; ++e) {
394          Value rem = (*_capacity)[e] - (*_flow)[e];
395          Node t = _graph.target(e);
396          if (_tolerance.positive(rem) && (*_pred)[t] == INVALID) {
397            _pred->set(t, e);
398            _queue[last++] = t;
399          }
400        }
401        for (InEdgeIt e(_graph, n); e != INVALID; ++e) {
402          Value rem = (*_flow)[e];
403          Node t = _graph.source(e);
404          if (_tolerance.positive(rem) && (*_pred)[t] == INVALID) {
405            _pred->set(t, e);
406            _queue[last++] = t;
407          }
408        }
409      }
410
411      if ((*_pred)[_target] != INVALID) {
412        Node n = _target;
413        Edge e = (*_pred)[n];
414
415        Value prem = (*_capacity)[e] - (*_flow)[e];
416        n = _graph.source(e);
417        while (n != _source) {
418          e = (*_pred)[n];
419          if (_graph.target(e) == n) {
420            Value rem = (*_capacity)[e] - (*_flow)[e];
421            if (rem < prem) prem = rem;
422            n = _graph.source(e);
423          } else {
424            Value rem = (*_flow)[e];
425            if (rem < prem) prem = rem;
426            n = _graph.target(e);   
427          }
428        }
429
430        n = _target;
431        e = (*_pred)[n];
432
433        _flow->set(e, (*_flow)[e] + prem);
434        n = _graph.source(e);
435        while (n != _source) {
436          e = (*_pred)[n];
437          if (_graph.target(e) == n) {
438            _flow->set(e, (*_flow)[e] + prem);
439            n = _graph.source(e);
440          } else {
441            _flow->set(e, (*_flow)[e] - prem);
442            n = _graph.target(e);   
443          }
444        }
445
446        _flow_value += prem;   
447        return true;
448      } else {
449        return false;
450      }
451    }
452
453    /// \brief Executes the algorithm
454    ///
455    /// It runs augmenting phases until the optimal solution is reached.
456    void start() {
457      while (augment()) {}
458    }
459
460    /// \brief runs the algorithm.
461    ///
462    /// It is just a shorthand for:
463    ///
464    ///\code
465    /// ek.init();
466    /// ek.start();
467    ///\endcode
468    void run() {
469      init();
470      start();
471    }
472
473    /// @}
474
475    /// \name Query Functions
476    /// The result of the Edmonds-Karp algorithm can be obtained using these
477    /// functions.\n
478    /// Before the use of these functions,
479    /// either run() or start() must be called.
480   
481    ///@{
482
483    /// \brief Returns the value of the maximum flow.
484    ///
485    /// Returns the value of the maximum flow by returning the excess
486    /// of the target node \c t. This value equals to the value of
487    /// the maximum flow already after the first phase.
488    Value flowValue() const {
489      return _flow_value;
490    }
491
492
493    /// \brief Returns the flow on the edge.
494    ///
495    /// Sets the \c flowMap to the flow on the edges. This method can
496    /// be called after the second phase of algorithm.
497    Value flow(const Edge& edge) const {
498      return (*_flow)[edge];
499    }
500
501    /// \brief Returns true when the node is on the source side of minimum cut.
502    ///
503
504    /// Returns true when the node is on the source side of minimum
505    /// cut. This method can be called both after running \ref
506    /// startFirstPhase() and \ref startSecondPhase().
507    bool minCut(const Node& node) const {
508      return (*_pred)[node] != INVALID;
509    }
510
511    /// \brief Returns a minimum value cut.
512    ///
513    /// Sets \c cut to the characteristic vector of a minimum value cut
514    /// It simply calls the minMinCut member.
515    /// \retval cut Write node bool map.
516    template <typename CutMap>
517    void minCutMap(CutMap& cutMap) const {
518      for (NodeIt n(_graph); n != INVALID; ++n) {
519        cutMap.set(n, (*_pred)[n] != INVALID);
520      }
521      cutMap.set(_source, true);
522    }   
523
524    /// @}
525
526  };
527
528}
529
530#endif
Note: See TracBrowser for help on using the repository browser.