lemon/edmonds_karp.h
author ladanyi
Fri, 14 Apr 2006 14:57:59 +0000
changeset 2048 b1a605b2f03c
parent 2036 9d0c8a205e58
child 2059 ebf3b2962554
permissions -rw-r--r--
make distcheck should work now
deba@2034
     1
/* -*- C++ -*-
deba@2034
     2
 *
deba@2034
     3
 * This file is a part of LEMON, a generic C++ optimization library
deba@2034
     4
 *
deba@2034
     5
 * Copyright (C) 2003-2006
deba@2034
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@2034
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@2034
     8
 *
deba@2034
     9
 * Permission to use, modify and distribute this software is granted
deba@2034
    10
 * provided that this copyright notice appears in all copies. For
deba@2034
    11
 * precise terms see the accompanying LICENSE file.
deba@2034
    12
 *
deba@2034
    13
 * This software is provided "AS IS" with no warranty of any kind,
deba@2034
    14
 * express or implied, and with no claim as to its suitability for any
deba@2034
    15
 * purpose.
deba@2034
    16
 *
deba@2034
    17
 */
deba@2034
    18
deba@2034
    19
#ifndef LEMON_EDMONDS_KARP_H
deba@2034
    20
#define LEMON_EDMONDS_KARP_H
deba@2034
    21
deba@2034
    22
/// \file
deba@2034
    23
/// \ingroup flowalgs
deba@2034
    24
/// \brief Implementation of the Edmonds-Karp algorithm.
deba@2034
    25
deba@2034
    26
#include <lemon/graph_adaptor.h>
deba@2034
    27
#include <lemon/tolerance.h>
deba@2034
    28
#include <lemon/bfs.h>
deba@2034
    29
deba@2034
    30
namespace lemon {
deba@2034
    31
deba@2034
    32
  /// \ingroup flowalgs
deba@2034
    33
  /// \brief Edmonds-Karp algorithms class.
deba@2034
    34
  ///
deba@2034
    35
  /// This class provides an implementation of the \e Edmonds-Karp \e
deba@2034
    36
  /// algorithm producing a flow of maximum value in a directed
deba@2034
    37
  /// graph. The Edmonds-Karp algorithm is slower than the Preflow algorithm
deba@2034
    38
  /// but it has an advantage of the step-by-step execution control with
deba@2034
    39
  /// feasible flow solutions. The \e source node, the \e target node, the \e
deba@2034
    40
  /// capacity of the edges and the \e starting \e flow value of the
deba@2034
    41
  /// edges should be passed to the algorithm through the
deba@2036
    42
  /// constructor.
deba@2034
    43
  ///
deba@2034
    44
  /// The time complexity of the algorithm is O(n * e^2) in worst case. 
deba@2034
    45
  /// Always try the preflow algorithm instead of this if you does not
deba@2034
    46
  /// have some additional reason than to compute the optimal flow which
deba@2034
    47
  /// has O(n^3) time complexity.
deba@2034
    48
  ///
deba@2034
    49
  /// \param _Graph The directed graph type the algorithm runs on.
deba@2034
    50
  /// \param _Number The number type of the capacities and the flow values.
deba@2034
    51
  /// \param _CapacityMap The capacity map type.
deba@2034
    52
  /// \param _FlowMap The flow map type.
deba@2034
    53
  /// \param _Tolerance The tolerance class to handle computation problems.
deba@2034
    54
  ///
deba@2034
    55
  /// \author Balazs Dezso 
deba@2034
    56
  template <typename _Graph, typename _Number,
deba@2034
    57
	    typename _CapacityMap = typename _Graph::template EdgeMap<_Number>,
deba@2034
    58
            typename _FlowMap = typename _Graph::template EdgeMap<_Number>,
deba@2034
    59
	    typename _Tolerance = Tolerance<_Number> >
deba@2034
    60
  class EdmondsKarp {
deba@2034
    61
  public:
deba@2034
    62
deba@2034
    63
    /// \brief \ref Exception for the case when the source equals the target.
deba@2034
    64
    ///
deba@2034
    65
    /// \ref Exception for the case when the source equals the target.
deba@2034
    66
    ///
deba@2034
    67
    class InvalidArgument : public lemon::LogicError {
deba@2034
    68
    public:
deba@2034
    69
      virtual const char* exceptionName() const {
deba@2034
    70
	return "lemon::EdmondsKarp::InvalidArgument";
deba@2034
    71
      }
deba@2034
    72
    };
deba@2034
    73
deba@2034
    74
deba@2034
    75
    /// \brief The graph type the algorithm runs on. 
deba@2034
    76
    typedef _Graph Graph;
deba@2034
    77
    /// \brief The value type of the algorithms.
deba@2034
    78
    typedef _Number Number;
deba@2034
    79
    /// \brief The capacity map on the edges.
deba@2034
    80
    typedef _CapacityMap CapacityMap;
deba@2034
    81
    /// \brief The flow map on the edges.
deba@2034
    82
    typedef _FlowMap FlowMap;
deba@2034
    83
    /// \brief The tolerance used by the algorithm.
deba@2034
    84
    typedef _Tolerance Tolerance;
deba@2034
    85
deba@2034
    86
    typedef ResGraphAdaptor<Graph, Number, CapacityMap, 
deba@2034
    87
                            FlowMap, Tolerance> ResGraph;
deba@2034
    88
deba@2034
    89
  private:
deba@2034
    90
deba@2034
    91
    typedef typename Graph::Node Node;
deba@2034
    92
    typedef typename Graph::Edge Edge;
deba@2034
    93
    
deba@2034
    94
    typedef typename Graph::NodeIt NodeIt;
deba@2034
    95
    typedef typename Graph::EdgeIt EdgeIt;
deba@2034
    96
    typedef typename Graph::InEdgeIt InEdgeIt;
deba@2034
    97
    typedef typename Graph::OutEdgeIt OutEdgeIt;
deba@2034
    98
    
deba@2034
    99
  public:
deba@2034
   100
deba@2034
   101
    /// \brief The constructor of the class.
deba@2034
   102
    ///
deba@2034
   103
    /// The constructor of the class. 
deba@2037
   104
    /// \param graph The directed graph the algorithm runs on. 
deba@2037
   105
    /// \param source The source node.
deba@2037
   106
    /// \param target The target node.
deba@2037
   107
    /// \param capacity The capacity of the edges. 
deba@2037
   108
    /// \param flow The flow of the edges. 
deba@2037
   109
    /// \param tolerance Tolerance class.
deba@2034
   110
    /// Except the graph, all of these parameters can be reset by
deba@2034
   111
    /// calling \ref source, \ref target, \ref capacityMap and \ref
deba@2034
   112
    /// flowMap, resp.
deba@2034
   113
    EdmondsKarp(const Graph& graph, Node source, Node target,
deba@2034
   114
                const CapacityMap& capacity, FlowMap& flow, 
deba@2034
   115
                const Tolerance& tolerance = Tolerance())
deba@2034
   116
      : _graph(graph), _capacity(capacity), _flow(flow), 
deba@2034
   117
        _tolerance(tolerance), _resgraph(graph, capacity, flow, tolerance), 
deba@2034
   118
        _source(source), _target(target) 
deba@2034
   119
    {
deba@2034
   120
      if (_source == _target) {
deba@2034
   121
        throw InvalidArgument();
deba@2034
   122
      }
deba@2034
   123
    }
deba@2034
   124
deba@2034
   125
    /// \brief Initializes the algorithm
deba@2034
   126
    /// 
deba@2034
   127
    /// It sets the flow to empty flow.
deba@2034
   128
    void init() {
deba@2034
   129
      for (EdgeIt it(_graph); it != INVALID; ++it) {
deba@2034
   130
        _flow.set(it, 0);
deba@2034
   131
      }
deba@2034
   132
      _value = 0;
deba@2034
   133
    }
deba@2034
   134
    
deba@2034
   135
    /// \brief Initializes the algorithm
deba@2034
   136
    /// 
deba@2034
   137
    /// If the flow map initially flow this let the flow map
deba@2034
   138
    /// unchanged but the flow value will be set by the flow
deba@2034
   139
    /// on the outedges from the source.
deba@2034
   140
    void flowInit() {
deba@2034
   141
      _value = 0;
deba@2034
   142
      for (OutEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
deba@2034
   143
        _value += _flow[jt];
deba@2034
   144
      }
deba@2034
   145
      for (InEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
deba@2034
   146
        _value -= _flow[jt];
deba@2034
   147
      }
deba@2034
   148
    }
deba@2034
   149
deba@2034
   150
    /// \brief Initializes the algorithm
deba@2034
   151
    /// 
deba@2034
   152
    /// If the flow map initially flow this let the flow map
deba@2034
   153
    /// unchanged but the flow value will be set by the flow
deba@2034
   154
    /// on the outedges from the source. It also checks that
deba@2034
   155
    /// the flow map really contains a flow.
deba@2034
   156
    /// \return %True when the flow map really a flow.
deba@2034
   157
    bool checkedFlowInit() {
deba@2034
   158
      _value = 0;
deba@2034
   159
      for (OutEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
deba@2034
   160
        _value += _flow[jt];
deba@2034
   161
      }
deba@2034
   162
      for (InEdgeIt jt(_graph, _source); jt != INVALID; ++jt) {
deba@2034
   163
        _value -= _flow[jt];
deba@2034
   164
      }
deba@2034
   165
      for (NodeIt it(_graph); it != INVALID; ++it) {
deba@2034
   166
        if (it == _source || it == _target) continue;
deba@2034
   167
        Number outFlow = 0;
deba@2034
   168
        for (OutEdgeIt jt(_graph, it); jt != INVALID; ++jt) {
deba@2034
   169
          outFlow += _flow[jt];
deba@2034
   170
        }
deba@2034
   171
        Number inFlow = 0;
deba@2034
   172
        for (InEdgeIt jt(_graph, it); jt != INVALID; ++jt) {
deba@2034
   173
          inFlow += _flow[jt];
deba@2034
   174
        }
deba@2034
   175
        if (_tolerance.different(outFlow, inFlow)) {
deba@2034
   176
          return false;
deba@2034
   177
        }
deba@2034
   178
      }
deba@2034
   179
      for (EdgeIt it(_graph); it != INVALID; ++it) {
deba@2034
   180
        if (_tolerance.less(_flow[it], 0)) return false;
deba@2034
   181
        if (_tolerance.less(_capacity[it], _flow[it])) return false;
deba@2034
   182
      }
deba@2034
   183
      return true;
deba@2034
   184
    }
deba@2034
   185
deba@2034
   186
    /// \brief Augment the solution on an edge shortest path.
deba@2034
   187
    /// 
deba@2034
   188
    /// Augment the solution on an edge shortest path. It search an
deba@2034
   189
    /// edge shortest path between the source and the target
deba@2034
   190
    /// in the residual graph with the bfs algoritm.
deba@2034
   191
    /// Then it increase the flow on this path with the minimal residual
deba@2034
   192
    /// capacity on the path. If there is not such path it gives back
deba@2034
   193
    /// false.
deba@2034
   194
    /// \return %False when the augmenting is not success so the
deba@2034
   195
    /// current flow is a feasible and optimal solution.
deba@2034
   196
    bool augment() {
deba@2034
   197
      typename Bfs<ResGraph>
deba@2034
   198
      ::template DefDistMap<NullMap<Node, int> >
deba@2034
   199
      ::Create bfs(_resgraph);
deba@2034
   200
deba@2034
   201
      NullMap<Node, int> distMap;
deba@2034
   202
      bfs.distMap(distMap);
deba@2034
   203
      
deba@2034
   204
      bfs.init();
deba@2034
   205
      bfs.addSource(_source);
deba@2034
   206
      bfs.start(_target);
deba@2034
   207
deba@2034
   208
      if (!bfs.reached(_target)) {
deba@2034
   209
        return false;
deba@2034
   210
      }
deba@2034
   211
      Number min = _resgraph.rescap(bfs.predEdge(_target));
deba@2034
   212
      for (Node it = bfs.predNode(_target); it != _source; 
deba@2034
   213
           it = bfs.predNode(it)) {
deba@2034
   214
        if (min > _resgraph.rescap(bfs.predEdge(it))) {
deba@2034
   215
          min = _resgraph.rescap(bfs.predEdge(it));
deba@2034
   216
        }
deba@2034
   217
      } 
deba@2034
   218
      for (Node it = _target; it != _source; it = bfs.predNode(it)) {
deba@2034
   219
        _resgraph.augment(bfs.predEdge(it), min);
deba@2034
   220
      }
deba@2034
   221
      _value += min;
deba@2034
   222
      return true;
deba@2034
   223
    }
deba@2034
   224
deba@2034
   225
    /// \brief Executes the algorithm
deba@2034
   226
    ///
deba@2034
   227
    /// It runs augmenting phases until the optimal solution is reached. 
deba@2034
   228
    void start() {
deba@2034
   229
      while (augment()) {}
deba@2034
   230
    }
deba@2034
   231
deba@2034
   232
    /// \brief Gives back the current flow value.
deba@2034
   233
    ///
deba@2034
   234
    /// Gives back the current flow _value.
deba@2034
   235
    Number flowValue() const {
deba@2034
   236
      return _value;
deba@2034
   237
    }
deba@2034
   238
deba@2034
   239
    /// \brief runs the algorithm.
deba@2034
   240
    /// 
deba@2034
   241
    /// It is just a shorthand for:
deba@2034
   242
    /// \code 
deba@2034
   243
    /// ek.init();
deba@2034
   244
    /// ek.start();
deba@2034
   245
    /// \endcode
deba@2034
   246
    void run() {
deba@2034
   247
      init();
deba@2034
   248
      start();
deba@2034
   249
    }
deba@2034
   250
deba@2034
   251
    /// \brief Returns a minimum value cut.
deba@2034
   252
    ///
deba@2034
   253
    /// Sets \c cut to the characteristic vector of a minimum value cut
deba@2034
   254
    /// It simply calls the minMinCut member.
deba@2037
   255
    /// \retval cut Write node bool map. 
deba@2034
   256
    template <typename CutMap>
deba@2034
   257
    void minCut(CutMap& cut) const {
deba@2034
   258
      minMinCut(cut);
deba@2034
   259
    }
deba@2034
   260
deba@2034
   261
    /// \brief Returns the inclusionwise minimum of the minimum value cuts.
deba@2034
   262
    ///
deba@2034
   263
    /// Sets \c cut to the characteristic vector of the minimum value cut
deba@2034
   264
    /// which is inclusionwise minimum. It is computed by processing a
deba@2034
   265
    /// bfs from the source node \c source in the residual graph.  
deba@2037
   266
    /// \retval cut Write node bool map. 
deba@2034
   267
    template <typename CutMap>
deba@2034
   268
    void minMinCut(CutMap& cut) const {
deba@2034
   269
deba@2034
   270
      typename Bfs<ResGraph>
deba@2034
   271
      ::template DefDistMap<NullMap<Node, int> >
deba@2034
   272
      ::template DefProcessedMap<CutMap>
deba@2034
   273
      ::Create bfs(_resgraph);
deba@2034
   274
deba@2034
   275
      NullMap<Node, int> distMap;
deba@2034
   276
      bfs.distMap(distMap);
deba@2034
   277
deba@2034
   278
      bfs.processedMap(cut);
deba@2034
   279
     
deba@2034
   280
      bfs.run(_source);
deba@2034
   281
    }
deba@2034
   282
deba@2034
   283
    /// \brief Returns the inclusionwise minimum of the minimum value cuts.
deba@2034
   284
    ///
deba@2034
   285
    /// Sets \c cut to the characteristic vector of the minimum value cut
deba@2034
   286
    /// which is inclusionwise minimum. It is computed by processing a
deba@2034
   287
    /// bfs from the source node \c source in the residual graph.  
deba@2037
   288
    /// \retval cut Write node bool map. 
deba@2034
   289
    template <typename CutMap>
deba@2034
   290
    void maxMinCut(CutMap& cut) const {
deba@2034
   291
deba@2034
   292
      typedef RevGraphAdaptor<const ResGraph> RevGraph;
deba@2034
   293
deba@2034
   294
      RevGraph revgraph(_resgraph);
deba@2034
   295
deba@2034
   296
      typename Bfs<RevGraph>
deba@2034
   297
      ::template DefDistMap<NullMap<Node, int> >
deba@2034
   298
      ::template DefPredMap<NullMap<Node, Edge> >
deba@2034
   299
      ::template DefProcessedMap<NotWriteMap<CutMap> >
deba@2034
   300
      ::Create bfs(revgraph);
deba@2034
   301
deba@2034
   302
      NullMap<Node, int> distMap;
deba@2034
   303
      bfs.distMap(distMap);
deba@2034
   304
deba@2034
   305
      NullMap<Node, Edge> predMap;
deba@2034
   306
      bfs.predMap(predMap);
deba@2034
   307
deba@2034
   308
      NotWriteMap<CutMap> notcut(cut);
deba@2034
   309
      bfs.processedMap(notcut);
deba@2034
   310
     
deba@2034
   311
      bfs.run(_target);
deba@2034
   312
    }
deba@2034
   313
deba@2034
   314
    /// \brief Returns the source node.
deba@2034
   315
    ///
deba@2034
   316
    /// Returns the source node.
deba@2034
   317
    /// 
deba@2034
   318
    Node source() const { 
deba@2034
   319
      return _source;
deba@2034
   320
    }
deba@2034
   321
deba@2034
   322
    /// \brief Returns the target node.
deba@2034
   323
    ///
deba@2034
   324
    /// Returns the target node.
deba@2034
   325
    /// 
deba@2034
   326
    Node target() const { 
deba@2034
   327
      return _target;
deba@2034
   328
    }
deba@2034
   329
deba@2034
   330
    /// \brief Returns a reference to capacity map.
deba@2034
   331
    ///
deba@2034
   332
    /// Returns a reference to capacity map.
deba@2034
   333
    /// 
deba@2034
   334
    const CapacityMap &capacityMap() const { 
deba@2034
   335
      return *_capacity;
deba@2034
   336
    }
deba@2034
   337
     
deba@2034
   338
    /// \brief Returns a reference to flow map.
deba@2034
   339
    ///
deba@2034
   340
    /// Returns a reference to flow map.
deba@2034
   341
    /// 
deba@2034
   342
    const FlowMap &flowMap() const { 
deba@2034
   343
      return *_flow;
deba@2034
   344
    }
deba@2034
   345
deba@2036
   346
    /// \brief Returns the tolerance used by algorithm.
deba@2036
   347
    ///
deba@2036
   348
    /// Returns the tolerance used by algorithm.
deba@2036
   349
    const Tolerance& tolerance() const {
deba@2036
   350
      return tolerance;
deba@2036
   351
    } 
deba@2036
   352
    
deba@2034
   353
  private:
deba@2034
   354
    
deba@2034
   355
    const Graph& _graph;
deba@2034
   356
    const CapacityMap& _capacity;
deba@2034
   357
    FlowMap& _flow;
deba@2034
   358
    Tolerance _tolerance;
deba@2034
   359
    
deba@2034
   360
    ResGraph _resgraph;
deba@2034
   361
    Node _source, _target;
deba@2034
   362
    Number _value;
deba@2034
   363
    
deba@2034
   364
  };
deba@2034
   365
deba@2034
   366
}
deba@2034
   367
deba@2034
   368
#endif