COIN-OR::LEMON - Graph Library

source: lemon/lemon/edmonds_karp.h @ 1224:92a884824429

Last change on this file since 1224:92a884824429 was 1224:92a884824429, checked in by Antal Nemes <thoneyvazul@…>, 9 years ago

Port Edmonds-Karp algorithm from svn -r3524 (#177)

File size: 13.9 KB
Line 
1/* -*- mode: C++; indent-tabs-mode: nil; -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library.
4 *
5 * Copyright (C) 2003-2010
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 GR Digraph type.
35  /// \param CAP Type of capacity map.
36  template <typename GR, typename CAP>
37  struct EdmondsKarpDefaultTraits {
38
39    /// \brief The digraph type the algorithm runs on.
40    typedef GR Digraph;
41
42    /// \brief The type of the map that stores the arc capacities.
43    ///
44    /// The type of the map that stores the arc capacities.
45    /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
46    typedef CAP CapacityMap;
47
48    /// \brief The type of the length of the arcs.
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 Digraph::template ArcMap<Value> FlowMap;
56
57    /// \brief Instantiates a FlowMap.
58    ///
59    /// This function instantiates a \ref FlowMap.
60    /// \param digraph The digraph, to which we would like to define the flow map.
61    static FlowMap* createFlowMap(const Digraph& digraph) {
62      return new FlowMap(digraph);
63    }
64
65    /// \brief The tolerance used by the algorithm
66    ///
67    /// The tolerance used by the algorithm to handle inexact computation.
68    typedef lemon::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 directed
78  /// digraphs. 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 arcs and the \e starting \e
82  /// flow value of the arcs 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 GR The digraph type the algorithm runs on.
90  /// \param CAP The capacity map type.
91  /// \param TR 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#ifdef DOXYGEN
97  template <typename GR, typename CAP, typename TR>
98#else
99  template <typename GR,
100            typename CAP = typename GR::template ArcMap<int>,
101            typename TR = EdmondsKarpDefaultTraits<GR, CAP> >
102#endif
103  class EdmondsKarp {
104  public:
105
106    typedef TR Traits;
107    typedef typename Traits::Digraph Digraph;
108    typedef typename Traits::CapacityMap CapacityMap;
109    typedef typename Traits::Value Value;
110
111    typedef typename Traits::FlowMap FlowMap;
112    typedef typename Traits::Tolerance Tolerance;
113
114  private:
115
116    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
117    typedef typename Digraph::template NodeMap<Arc> PredMap;
118   
119    const Digraph& _graph;
120    const CapacityMap* _capacity;
121
122    Node _source, _target;
123
124    FlowMap* _flow;
125    bool _local_flow;
126
127    PredMap* _pred;
128    std::vector<Node> _queue;
129   
130    Tolerance _tolerance;
131    Value _flow_value;
132
133    void createStructures() {
134      if (!_flow) {
135        _flow = Traits::createFlowMap(_graph);
136        _local_flow = true;
137      }
138      if (!_pred) {
139        _pred = new PredMap(_graph);
140      }
141      _queue.resize(countNodes(_graph));
142    }
143
144    void destroyStructures() {
145      if (_local_flow) {
146        delete _flow;
147      }
148      if (_pred) {
149        delete _pred;
150      }
151    }
152   
153  public:
154
155    ///\name Named template parameters
156
157    ///@{
158
159    template <typename T>
160    struct DefFlowMapTraits : public Traits {
161      typedef T FlowMap;
162      static FlowMap *createFlowMap(const Digraph&) {
163        LEMON_ASSERT(false,"Uninitialized parameter.");
164        return 0;
165      }
166    };
167
168    /// \brief \ref named-templ-param "Named parameter" for setting
169    /// FlowMap type
170    ///
171    /// \ref named-templ-param "Named parameter" for setting FlowMap
172    /// type
173    template <typename T>
174    struct DefFlowMap
175      : public EdmondsKarp<Digraph, CapacityMap, DefFlowMapTraits<T> > {
176      typedef EdmondsKarp<Digraph, CapacityMap, DefFlowMapTraits<T> >
177      Create;
178    };
179
180
181    /// @}
182
183  protected:
184   
185    EdmondsKarp() {}
186
187  public:
188
189    /// \brief The constructor of the class.
190    ///
191    /// The constructor of the class.
192    /// \param digraph The digraph the algorithm runs on.
193    /// \param capacity The capacity of the arcs.
194    /// \param source The source node.
195    /// \param target The target node.
196    EdmondsKarp(const Digraph& digraph, const CapacityMap& capacity,
197                Node source, Node target)
198      : _graph(digraph), _capacity(&capacity), _source(source), _target(target),
199        _flow(0), _local_flow(false), _pred(0), _tolerance(), _flow_value()
200    {
201      LEMON_ASSERT(_source != _target,"Flow source and target are the same nodes.");
202    }
203
204    /// \brief Destructor.
205    ///
206    /// Destructor.
207    ~EdmondsKarp() {
208      destroyStructures();
209    }
210
211    /// \brief Sets the capacity map.
212    ///
213    /// Sets the capacity map.
214    /// \return \c (*this)
215    EdmondsKarp& capacityMap(const CapacityMap& map) {
216      _capacity = &map;
217      return *this;
218    }
219
220    /// \brief Sets the flow map.
221    ///
222    /// Sets the flow map.
223    /// \return \c (*this)
224    EdmondsKarp& flowMap(FlowMap& map) {
225      if (_local_flow) {
226        delete _flow;
227        _local_flow = false;
228      }
229      _flow = &map;
230      return *this;
231    }
232
233    /// \brief Returns the flow map.
234    ///
235    /// \return The flow map.
236    const FlowMap& flowMap() const {
237      return *_flow;
238    }
239
240    /// \brief Sets the source node.
241    ///
242    /// Sets the source node.
243    /// \return \c (*this)
244    EdmondsKarp& source(const Node& node) {
245      _source = node;
246      return *this;
247    }
248
249    /// \brief Sets the target node.
250    ///
251    /// Sets the target node.
252    /// \return \c (*this)
253    EdmondsKarp& target(const Node& node) {
254      _target = node;
255      return *this;
256    }
257
258    /// \brief Sets the tolerance used by algorithm.
259    ///
260    /// Sets the tolerance used by algorithm.
261    EdmondsKarp& tolerance(const Tolerance& tolerance) {
262      _tolerance = tolerance;
263      return *this;
264    }
265
266    /// \brief Returns the tolerance used by algorithm.
267    ///
268    /// Returns the tolerance used by algorithm.
269    const Tolerance& tolerance() const {
270      return _tolerance;
271    }
272
273    /// \name Execution control
274    /// The simplest way to execute the
275    /// algorithm is to use the \c run() member functions.
276    /// \n
277    /// If you need more control on initial solution or
278    /// execution then you have to call one \ref init() function and then
279    /// the start() or multiple times the \c augment() member function. 
280   
281    ///@{
282
283    /// \brief Initializes the algorithm
284    ///
285    /// Sets the flow to empty flow.
286    void init() {
287      createStructures();
288      for (ArcIt it(_graph); it != INVALID; ++it) {
289        _flow->set(it, 0);
290      }
291      _flow_value = 0;
292    }
293   
294    /// \brief Initializes the algorithm
295    ///
296    /// Initializes the flow to the \c flowMap. The \c flowMap should
297    /// contain a feasible flow, ie. in each node excluding the source
298    /// and the target the incoming flow should be equal to the
299    /// outgoing flow.
300    template <typename FlowMap>
301    void flowInit(const FlowMap& flowMap) {
302      createStructures();
303      for (ArcIt e(_graph); e != INVALID; ++e) {
304        _flow->set(e, flowMap[e]);
305      }
306      _flow_value = 0;
307      for (OutArcIt jt(_graph, _source); jt != INVALID; ++jt) {
308        _flow_value += (*_flow)[jt];
309      }
310      for (InArcIt jt(_graph, _source); jt != INVALID; ++jt) {
311        _flow_value -= (*_flow)[jt];
312      }
313    }
314
315    /// \brief Initializes the algorithm
316    ///
317    /// Initializes the flow to the \c flowMap. The \c flowMap should
318    /// contain a feasible flow, ie. in each node excluding the source
319    /// and the target the incoming flow should be equal to the
320    /// outgoing flow. 
321    /// \return %False when the given flowMap does not contain
322    /// feasible flow.
323    template <typename FlowMap>
324    bool checkedFlowInit(const FlowMap& flowMap) {
325      createStructures();
326      for (ArcIt e(_graph); e != INVALID; ++e) {
327        _flow->set(e, flowMap[e]);
328      }
329      for (NodeIt it(_graph); it != INVALID; ++it) {
330        if (it == _source || it == _target) continue;
331        Value outFlow = 0;
332        for (OutArcIt jt(_graph, it); jt != INVALID; ++jt) {
333          outFlow += (*_flow)[jt];
334        }
335        Value inFlow = 0;
336        for (InArcIt jt(_graph, it); jt != INVALID; ++jt) {
337          inFlow += (*_flow)[jt];
338        }
339        if (_tolerance.different(outFlow, inFlow)) {
340          return false;
341        }
342      }
343      for (ArcIt it(_graph); it != INVALID; ++it) {
344        if (_tolerance.less((*_flow)[it], 0)) return false;
345        if (_tolerance.less((*_capacity)[it], (*_flow)[it])) return false;
346      }
347      _flow_value = 0;
348      for (OutArcIt jt(_graph, _source); jt != INVALID; ++jt) {
349        _flow_value += (*_flow)[jt];
350      }
351      for (InArcIt jt(_graph, _source); jt != INVALID; ++jt) {
352        _flow_value -= (*_flow)[jt];
353      }
354      return true;
355    }
356
357    /// \brief Augment the solution on an arc shortest path.
358    ///
359    /// Augment the solution on an arc shortest path. It searches an
360    /// arc shortest path between the source and the target
361    /// in the residual digraph by the bfs algoritm.
362    /// Then it increases the flow on this path with the minimal residual
363    /// capacity on the path. If there is no such path it gives back
364    /// false.
365    /// \return %False when the augmenting didn't success so the
366    /// current flow is a feasible and optimal solution.
367    bool augment() {
368      for (NodeIt n(_graph); n != INVALID; ++n) {
369        _pred->set(n, INVALID);
370      }
371     
372      int first = 0, last = 1;
373     
374      _queue[0] = _source;
375      _pred->set(_source, OutArcIt(_graph, _source));
376
377      while (first != last && (*_pred)[_target] == INVALID) {
378        Node n = _queue[first++];
379       
380        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
381          Value rem = (*_capacity)[e] - (*_flow)[e];
382          Node t = _graph.target(e);
383          if (_tolerance.positive(rem) && (*_pred)[t] == INVALID) {
384            _pred->set(t, e);
385            _queue[last++] = t;
386          }
387        }
388        for (InArcIt e(_graph, n); e != INVALID; ++e) {
389          Value rem = (*_flow)[e];
390          Node t = _graph.source(e);
391          if (_tolerance.positive(rem) && (*_pred)[t] == INVALID) {
392            _pred->set(t, e);
393            _queue[last++] = t;
394          }
395        }
396      }
397
398      if ((*_pred)[_target] != INVALID) {
399        Node n = _target;
400        Arc e = (*_pred)[n];
401
402        Value prem = (*_capacity)[e] - (*_flow)[e];
403        n = _graph.source(e);
404        while (n != _source) {
405          e = (*_pred)[n];
406          if (_graph.target(e) == n) {
407            Value rem = (*_capacity)[e] - (*_flow)[e];
408            if (rem < prem) prem = rem;
409            n = _graph.source(e);
410          } else {
411            Value rem = (*_flow)[e];
412            if (rem < prem) prem = rem;
413            n = _graph.target(e);   
414          }
415        }
416
417        n = _target;
418        e = (*_pred)[n];
419
420        _flow->set(e, (*_flow)[e] + prem);
421        n = _graph.source(e);
422        while (n != _source) {
423          e = (*_pred)[n];
424          if (_graph.target(e) == n) {
425            _flow->set(e, (*_flow)[e] + prem);
426            n = _graph.source(e);
427          } else {
428            _flow->set(e, (*_flow)[e] - prem);
429            n = _graph.target(e);   
430          }
431        }
432
433        _flow_value += prem;   
434        return true;
435      } else {
436        return false;
437      }
438    }
439
440    /// \brief Executes the algorithm
441    ///
442    /// It runs augmenting phases until the optimal solution is reached.
443    void start() {
444      while (augment()) {}
445    }
446
447    /// \brief Runs the algorithm.
448    ///
449    /// It is just a shorthand for:
450    ///
451    ///\code
452    /// ek.init();
453    /// ek.start();
454    ///\endcode
455    void run() {
456      init();
457      start();
458    }
459
460    /// @}
461
462    /// \name Query Functions
463    /// The result of the Edmonds-Karp algorithm can be obtained using these
464    /// functions.\n
465    /// Before the use of these functions,
466    /// either run() or start() must be called.
467   
468    ///@{
469
470    /// \brief Returns the value of the maximum flow.
471    ///
472    /// Returns the value of the maximum flow by returning the excess
473    /// of the target node \c t.
474
475    Value flowValue() const {
476      return _flow_value;
477    }
478
479
480    /// \brief Returns the flow on the arc.
481    ///
482    /// Sets the \c flowMap to the flow on the arcs.
483    Value flow(const Arc& arc) const {
484      return (*_flow)[arc];
485    }
486
487    /// \brief Returns true when the node is on the source side of minimum cut.
488    ///
489
490    /// Returns true when the node is on the source side of minimum
491    /// cut.
492
493    bool minCut(const Node& node) const {
494      return ((*_pred)[node] != INVALID) or node == _source;
495    }
496
497    /// \brief Returns a minimum value cut.
498    ///
499    /// Sets \c cutMap to the characteristic vector of a minimum value cut.
500
501    template <typename CutMap>
502    void minCutMap(CutMap& cutMap) const {
503      for (NodeIt n(_graph); n != INVALID; ++n) {
504        cutMap.set(n, (*_pred)[n] != INVALID);
505      }
506      cutMap.set(_source, true);
507    }   
508
509    /// @}
510
511  };
512
513}
514
515#endif
Note: See TracBrowser for help on using the repository browser.