deba@2034: /* -*- C++ -*- deba@2034: * deba@2034: * This file is a part of LEMON, a generic C++ optimization library deba@2034: * alpar@2553: * Copyright (C) 2003-2008 deba@2034: * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport deba@2034: * (Egervary Research Group on Combinatorial Optimization, EGRES). deba@2034: * deba@2034: * Permission to use, modify and distribute this software is granted deba@2034: * provided that this copyright notice appears in all copies. For deba@2034: * precise terms see the accompanying LICENSE file. deba@2034: * deba@2034: * This software is provided "AS IS" with no warranty of any kind, deba@2034: * express or implied, and with no claim as to its suitability for any deba@2034: * purpose. deba@2034: * deba@2034: */ deba@2034: deba@2034: #ifndef LEMON_EDMONDS_KARP_H deba@2034: #define LEMON_EDMONDS_KARP_H deba@2034: deba@2034: /// \file deba@2376: /// \ingroup max_flow deba@2034: /// \brief Implementation of the Edmonds-Karp algorithm. deba@2034: deba@2034: #include deba@2514: #include deba@2034: deba@2034: namespace lemon { deba@2034: deba@2514: /// \brief Default traits class of EdmondsKarp class. deba@2514: /// deba@2514: /// Default traits class of EdmondsKarp class. deba@2514: /// \param _Graph Graph type. deba@2514: /// \param _CapacityMap Type of capacity map. deba@2514: template deba@2514: struct EdmondsKarpDefaultTraits { deba@2514: deba@2514: /// \brief The graph type the algorithm runs on. deba@2514: typedef _Graph Graph; deba@2514: deba@2514: /// \brief The type of the map that stores the edge capacities. deba@2514: /// deba@2514: /// The type of the map that stores the edge capacities. deba@2514: /// It must meet the \ref concepts::ReadMap "ReadMap" concept. deba@2514: typedef _CapacityMap CapacityMap; deba@2514: deba@2514: /// \brief The type of the length of the edges. deba@2514: typedef typename CapacityMap::Value Value; deba@2514: deba@2514: /// \brief The map type that stores the flow values. deba@2514: /// deba@2514: /// The map type that stores the flow values. deba@2514: /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept. deba@2514: typedef typename Graph::template EdgeMap FlowMap; deba@2514: deba@2514: /// \brief Instantiates a FlowMap. deba@2514: /// deba@2514: /// This function instantiates a \ref FlowMap. deba@2514: /// \param graph The graph, to which we would like to define the flow map. deba@2514: static FlowMap* createFlowMap(const Graph& graph) { deba@2514: return new FlowMap(graph); deba@2514: } deba@2514: deba@2514: /// \brief The tolerance used by the algorithm deba@2514: /// deba@2514: /// The tolerance used by the algorithm to handle inexact computation. deba@2514: typedef Tolerance Tolerance; deba@2514: deba@2514: }; deba@2514: deba@2376: /// \ingroup max_flow deba@2514: /// deba@2034: /// \brief Edmonds-Karp algorithms class. deba@2034: /// deba@2034: /// This class provides an implementation of the \e Edmonds-Karp \e deba@2034: /// algorithm producing a flow of maximum value in a directed deba@2514: /// graphs. The Edmonds-Karp algorithm is slower than the Preflow deba@2514: /// algorithm but it has an advantage of the step-by-step execution deba@2514: /// control with feasible flow solutions. The \e source node, the \e deba@2514: /// target node, the \e capacity of the edges and the \e starting \e deba@2514: /// flow value of the edges should be passed to the algorithm deba@2514: /// through the constructor. deba@2034: /// deba@2514: /// The time complexity of the algorithm is \f$ O(nm^2) \f$ in deba@2059: /// worst case. Always try the preflow algorithm instead of this if deba@2113: /// you just want to compute the optimal flow. deba@2034: /// deba@2034: /// \param _Graph The directed graph type the algorithm runs on. deba@2034: /// \param _CapacityMap The capacity map type. deba@2514: /// \param _Traits Traits class to set various data types used by deba@2514: /// the algorithm. The default traits class is \ref deba@2514: /// EdmondsKarpDefaultTraits. See \ref EdmondsKarpDefaultTraits for the deba@2514: /// documentation of a Edmonds-Karp traits class. deba@2034: /// deba@2034: /// \author Balazs Dezso deba@2059: #ifdef DOXYGEN deba@2514: template deba@2514: #else deba@2514: template , deba@2514: typename _Traits = EdmondsKarpDefaultTraits<_Graph, _CapacityMap> > deba@2059: #endif deba@2034: class EdmondsKarp { deba@2034: public: deba@2034: deba@2514: typedef _Traits Traits; deba@2514: typedef typename Traits::Graph Graph; deba@2514: typedef typename Traits::CapacityMap CapacityMap; deba@2514: typedef typename Traits::Value Value; deba@2514: deba@2514: typedef typename Traits::FlowMap FlowMap; deba@2514: typedef typename Traits::Tolerance Tolerance; deba@2514: deba@2034: /// \brief \ref Exception for the case when the source equals the target. deba@2034: /// deba@2034: /// \ref Exception for the case when the source equals the target. deba@2034: /// deba@2034: class InvalidArgument : public lemon::LogicError { deba@2034: public: alpar@2151: virtual const char* what() const throw() { deba@2034: return "lemon::EdmondsKarp::InvalidArgument"; deba@2034: } deba@2034: }; deba@2034: deba@2034: deba@2034: private: deba@2034: deba@2514: GRAPH_TYPEDEFS(typename Graph); deba@2514: typedef typename Graph::template NodeMap PredMap; deba@2034: deba@2514: const Graph& _graph; deba@2514: const CapacityMap* _capacity; deba@2514: deba@2514: Node _source, _target; deba@2514: deba@2514: FlowMap* _flow; deba@2514: bool _local_flow; deba@2514: deba@2514: PredMap* _pred; deba@2514: std::vector _queue; deba@2514: deba@2514: Tolerance _tolerance; deba@2514: Value _flow_value; deba@2514: deba@2514: void createStructures() { deba@2514: if (!_flow) { deba@2514: _flow = Traits::createFlowMap(_graph); deba@2514: _local_flow = true; deba@2514: } deba@2514: if (!_pred) { deba@2514: _pred = new PredMap(_graph); deba@2514: } deba@2514: _queue.resize(countNodes(_graph)); deba@2514: } deba@2514: deba@2514: void destroyStructures() { deba@2514: if (_local_flow) { deba@2514: delete _flow; deba@2514: } deba@2514: if (_pred) { deba@2514: delete _pred; deba@2514: } deba@2514: } deba@2034: deba@2034: public: deba@2034: deba@2514: ///\name Named template parameters deba@2514: deba@2514: ///@{ deba@2514: deba@2514: template deba@2514: struct DefFlowMapTraits : public Traits { deba@2514: typedef _FlowMap FlowMap; deba@2514: static FlowMap *createFlowMap(const Graph&) { deba@2514: throw UninitializedParameter(); deba@2514: } deba@2514: }; deba@2514: deba@2514: /// \brief \ref named-templ-param "Named parameter" for setting deba@2514: /// FlowMap type deba@2514: /// deba@2514: /// \ref named-templ-param "Named parameter" for setting FlowMap deba@2514: /// type deba@2514: template deba@2514: struct DefFlowMap deba@2514: : public EdmondsKarp > { deba@2514: typedef EdmondsKarp > deba@2514: Create; deba@2514: }; deba@2514: deba@2514: deba@2514: /// @} deba@2514: deba@2527: protected: deba@2527: deba@2527: EdmondsKarp() {} deba@2527: deba@2527: public: deba@2527: deba@2034: /// \brief The constructor of the class. deba@2034: /// deba@2034: /// The constructor of the class. deba@2037: /// \param graph The directed graph the algorithm runs on. deba@2514: /// \param capacity The capacity of the edges. deba@2037: /// \param source The source node. deba@2037: /// \param target The target node. deba@2514: EdmondsKarp(const Graph& graph, const CapacityMap& capacity, deba@2514: Node source, Node target) deba@2514: : _graph(graph), _capacity(&capacity), _source(source), _target(target), deba@2514: _flow(0), _local_flow(false), _pred(0), _tolerance(), _flow_value() deba@2034: { deba@2034: if (_source == _target) { deba@2034: throw InvalidArgument(); deba@2034: } deba@2034: } deba@2034: deba@2514: /// \brief Destrcutor. deba@2514: /// deba@2514: /// Destructor. deba@2514: ~EdmondsKarp() { deba@2514: destroyStructures(); deba@2514: } deba@2514: deba@2514: /// \brief Sets the capacity map. deba@2514: /// deba@2514: /// Sets the capacity map. deba@2514: /// \return \c (*this) deba@2514: EdmondsKarp& capacityMap(const CapacityMap& map) { deba@2514: _capacity = ↦ deba@2514: return *this; deba@2514: } deba@2514: deba@2514: /// \brief Sets the flow map. deba@2514: /// deba@2514: /// Sets the flow map. deba@2514: /// \return \c (*this) deba@2514: EdmondsKarp& flowMap(FlowMap& map) { deba@2514: if (_local_flow) { deba@2514: delete _flow; deba@2514: _local_flow = false; deba@2514: } deba@2514: _flow = ↦ deba@2514: return *this; deba@2514: } deba@2514: deba@2514: /// \brief Returns the flow map. deba@2514: /// deba@2514: /// \return The flow map. deba@2514: const FlowMap& flowMap() { deba@2514: return *_flow; deba@2514: } deba@2514: deba@2514: /// \brief Sets the source node. deba@2514: /// deba@2514: /// Sets the source node. deba@2514: /// \return \c (*this) deba@2514: EdmondsKarp& source(const Node& node) { deba@2514: _source = node; deba@2514: return *this; deba@2514: } deba@2514: deba@2514: /// \brief Sets the target node. deba@2514: /// deba@2514: /// Sets the target node. deba@2514: /// \return \c (*this) deba@2514: EdmondsKarp& target(const Node& node) { deba@2514: _target = node; deba@2514: return *this; deba@2514: } deba@2514: deba@2514: /// \brief Sets the tolerance used by algorithm. deba@2514: /// deba@2514: /// Sets the tolerance used by algorithm. deba@2514: EdmondsKarp& tolerance(const Tolerance& tolerance) const { deba@2514: _tolerance = tolerance; deba@2514: return *this; deba@2514: } deba@2514: deba@2514: /// \brief Returns the tolerance used by algorithm. deba@2514: /// deba@2514: /// Returns the tolerance used by algorithm. deba@2514: const Tolerance& tolerance() const { deba@2514: return tolerance; deba@2514: } deba@2514: deba@2514: /// \name Execution control The simplest way to execute the deba@2514: /// algorithm is to use the \c run() member functions. deba@2514: /// \n deba@2514: /// If you need more control on initial solution or deba@2514: /// execution then you have to call one \ref init() function and then deba@2514: /// the start() or multiple times the \c augment() member function. deba@2514: deba@2514: ///@{ deba@2514: deba@2034: /// \brief Initializes the algorithm deba@2034: /// deba@2034: /// It sets the flow to empty flow. deba@2034: void init() { deba@2514: createStructures(); deba@2034: for (EdgeIt it(_graph); it != INVALID; ++it) { deba@2514: _flow->set(it, 0); deba@2034: } deba@2514: _flow_value = 0; deba@2034: } deba@2034: deba@2034: /// \brief Initializes the algorithm deba@2034: /// deba@2514: /// Initializes the flow to the \c flowMap. The \c flowMap should deba@2514: /// contain a feasible flow, ie. in each node excluding the source deba@2514: /// and the target the incoming flow should be equal to the deba@2514: /// outgoing flow. deba@2514: template deba@2514: void flowInit(const FlowMap& flowMap) { deba@2514: createStructures(); deba@2514: for (EdgeIt e(_graph); e != INVALID; ++e) { deba@2514: _flow->set(e, flowMap[e]); deba@2514: } deba@2514: _flow_value = 0; deba@2034: for (OutEdgeIt jt(_graph, _source); jt != INVALID; ++jt) { deba@2514: _flow_value += (*_flow)[jt]; deba@2034: } deba@2034: for (InEdgeIt jt(_graph, _source); jt != INVALID; ++jt) { deba@2514: _flow_value -= (*_flow)[jt]; deba@2034: } deba@2034: } deba@2034: deba@2034: /// \brief Initializes the algorithm deba@2034: /// deba@2514: /// Initializes the flow to the \c flowMap. The \c flowMap should deba@2514: /// contain a feasible flow, ie. in each node excluding the source deba@2514: /// and the target the incoming flow should be equal to the deba@2514: /// outgoing flow. deba@2514: /// \return %False when the given flowMap does not contain deba@2514: /// feasible flow. deba@2514: template deba@2514: bool checkedFlowInit(const FlowMap& flowMap) { deba@2514: createStructures(); deba@2514: for (EdgeIt e(_graph); e != INVALID; ++e) { deba@2514: _flow->set(e, flowMap[e]); deba@2034: } deba@2034: for (NodeIt it(_graph); it != INVALID; ++it) { deba@2034: if (it == _source || it == _target) continue; deba@2514: Value outFlow = 0; deba@2034: for (OutEdgeIt jt(_graph, it); jt != INVALID; ++jt) { deba@2514: outFlow += (*_flow)[jt]; deba@2034: } deba@2514: Value inFlow = 0; deba@2034: for (InEdgeIt jt(_graph, it); jt != INVALID; ++jt) { deba@2514: inFlow += (*_flow)[jt]; deba@2034: } deba@2034: if (_tolerance.different(outFlow, inFlow)) { deba@2034: return false; deba@2034: } deba@2034: } deba@2034: for (EdgeIt it(_graph); it != INVALID; ++it) { deba@2514: if (_tolerance.less((*_flow)[it], 0)) return false; deba@2514: if (_tolerance.less((*_capacity)[it], (*_flow)[it])) return false; deba@2514: } deba@2514: _flow_value = 0; deba@2514: for (OutEdgeIt jt(_graph, _source); jt != INVALID; ++jt) { deba@2514: _flow_value += (*_flow)[jt]; deba@2514: } deba@2514: for (InEdgeIt jt(_graph, _source); jt != INVALID; ++jt) { deba@2514: _flow_value -= (*_flow)[jt]; deba@2034: } deba@2034: return true; deba@2034: } deba@2034: deba@2034: /// \brief Augment the solution on an edge shortest path. deba@2034: /// deba@2034: /// Augment the solution on an edge shortest path. It search an deba@2034: /// edge shortest path between the source and the target deba@2034: /// in the residual graph with the bfs algoritm. deba@2034: /// Then it increase the flow on this path with the minimal residual deba@2034: /// capacity on the path. If there is not such path it gives back deba@2034: /// false. deba@2034: /// \return %False when the augmenting is not success so the deba@2034: /// current flow is a feasible and optimal solution. deba@2034: bool augment() { deba@2514: for (NodeIt n(_graph); n != INVALID; ++n) { deba@2514: _pred->set(n, INVALID); deba@2514: } deba@2514: deba@2514: int first = 0, last = 1; deba@2514: deba@2514: _queue[0] = _source; deba@2514: _pred->set(_source, OutEdgeIt(_graph, _source)); deba@2034: deba@2514: while (first != last && (*_pred)[_target] == INVALID) { deba@2514: Node n = _queue[first++]; deba@2514: deba@2514: for (OutEdgeIt e(_graph, n); e != INVALID; ++e) { deba@2514: Value rem = (*_capacity)[e] - (*_flow)[e]; deba@2514: Node t = _graph.target(e); deba@2514: if (_tolerance.positive(rem) && (*_pred)[t] == INVALID) { deba@2514: _pred->set(t, e); deba@2514: _queue[last++] = t; deba@2514: } deba@2514: } deba@2514: for (InEdgeIt e(_graph, n); e != INVALID; ++e) { deba@2514: Value rem = (*_flow)[e]; deba@2514: Node t = _graph.source(e); deba@2514: if (_tolerance.positive(rem) && (*_pred)[t] == INVALID) { deba@2514: _pred->set(t, e); deba@2514: _queue[last++] = t; deba@2514: } deba@2514: } deba@2514: } deba@2034: deba@2514: if ((*_pred)[_target] != INVALID) { deba@2514: Node n = _target; deba@2514: Edge e = (*_pred)[n]; deba@2514: deba@2514: Value prem = (*_capacity)[e] - (*_flow)[e]; deba@2514: n = _graph.source(e); deba@2514: while (n != _source) { deba@2514: e = (*_pred)[n]; deba@2514: if (_graph.target(e) == n) { deba@2514: Value rem = (*_capacity)[e] - (*_flow)[e]; deba@2514: if (rem < prem) prem = rem; deba@2514: n = _graph.source(e); deba@2514: } else { deba@2514: Value rem = (*_flow)[e]; deba@2514: if (rem < prem) prem = rem; deba@2514: n = _graph.target(e); deba@2514: } deba@2514: } deba@2514: deba@2514: n = _target; deba@2514: e = (*_pred)[n]; deba@2514: deba@2514: _flow->set(e, (*_flow)[e] + prem); deba@2514: n = _graph.source(e); deba@2514: while (n != _source) { deba@2514: e = (*_pred)[n]; deba@2514: if (_graph.target(e) == n) { deba@2514: _flow->set(e, (*_flow)[e] + prem); deba@2514: n = _graph.source(e); deba@2514: } else { deba@2514: _flow->set(e, (*_flow)[e] - prem); deba@2514: n = _graph.target(e); deba@2514: } deba@2514: } deba@2514: deba@2514: _flow_value += prem; deba@2514: return true; deba@2514: } else { deba@2514: return false; deba@2034: } deba@2034: } deba@2034: deba@2034: /// \brief Executes the algorithm deba@2034: /// deba@2034: /// It runs augmenting phases until the optimal solution is reached. deba@2034: void start() { deba@2034: while (augment()) {} deba@2034: } deba@2034: deba@2034: /// \brief runs the algorithm. deba@2034: /// deba@2034: /// It is just a shorthand for: deba@2059: /// deba@2059: ///\code deba@2034: /// ek.init(); deba@2034: /// ek.start(); deba@2059: ///\endcode deba@2034: void run() { deba@2034: init(); deba@2034: start(); deba@2034: } deba@2034: deba@2514: /// @} deba@2514: deba@2514: /// \name Query Functions deba@2522: /// The result of the Edmonds-Karp algorithm can be obtained using these deba@2514: /// functions.\n deba@2514: /// Before the use of these functions, deba@2514: /// either run() or start() must be called. deba@2514: deba@2514: ///@{ deba@2514: deba@2514: /// \brief Returns the value of the maximum flow. deba@2514: /// deba@2514: /// Returns the value of the maximum flow by returning the excess deba@2514: /// of the target node \c t. This value equals to the value of deba@2514: /// the maximum flow already after the first phase. deba@2514: Value flowValue() const { deba@2514: return _flow_value; deba@2514: } deba@2514: deba@2514: deba@2514: /// \brief Returns the flow on the edge. deba@2514: /// deba@2514: /// Sets the \c flowMap to the flow on the edges. This method can deba@2514: /// be called after the second phase of algorithm. deba@2514: Value flow(const Edge& edge) const { deba@2514: return (*_flow)[edge]; deba@2514: } deba@2514: deba@2514: /// \brief Returns true when the node is on the source side of minimum cut. deba@2514: /// deba@2514: deba@2514: /// Returns true when the node is on the source side of minimum deba@2514: /// cut. This method can be called both after running \ref deba@2514: /// startFirstPhase() and \ref startSecondPhase(). deba@2514: bool minCut(const Node& node) const { deba@2514: return (*_pred)[node] != INVALID; deba@2514: } deba@2514: deba@2034: /// \brief Returns a minimum value cut. deba@2034: /// deba@2034: /// Sets \c cut to the characteristic vector of a minimum value cut deba@2034: /// It simply calls the minMinCut member. deba@2037: /// \retval cut Write node bool map. deba@2034: template deba@2514: void minCutMap(CutMap& cutMap) const { deba@2514: for (NodeIt n(_graph); n != INVALID; ++n) { deba@2514: cutMap.set(n, (*_pred)[n] != INVALID); deba@2514: } deba@2514: cutMap.set(_source, true); deba@2514: } deba@2034: deba@2514: /// @} deba@2034: deba@2034: }; deba@2034: deba@2034: } deba@2034: deba@2034: #endif