Some bug fix.
Added: Concept check for maps and IdReader/Writer classes
Some runtime check
2 * lemon/min_cost_flow.h - Part of LEMON, a generic C++ optimization library
4 * Copyright (C) 2005 Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
5 * (Egervary Research Group on Combinatorial Optimization, EGRES).
7 * Permission to use, modify and distribute this software is granted
8 * provided that this copyright notice appears in all copies. For
9 * precise terms see the accompanying LICENSE file.
11 * This software is provided "AS IS" with no warranty of any kind,
12 * express or implied, and with no claim as to its suitability for any
17 #ifndef LEMON_MIN_COST_FLOW_H
18 #define LEMON_MIN_COST_FLOW_H
22 ///\brief An algorithm for finding a flow of value \c k (for small values of \c k) having minimal total cost
25 #include <lemon/dijkstra.h>
26 #include <lemon/graph_adaptor.h>
27 #include <lemon/maps.h>
32 /// \addtogroup flowalgs
35 ///\brief Implementation of an algorithm for finding a flow of value \c k
36 ///(for small values of \c k) having minimal total cost between 2 nodes
39 /// The class \ref lemon::MinCostFlow "MinCostFlow" implements an
40 /// algorithm for finding a flow of value \c k having minimal total
41 /// cost from a given source node to a given target node in an
42 /// edge-weighted directed graph. To this end, the edge-capacities
43 /// and edge-weights have to be nonnegative. The edge-capacities
44 /// should be integers, but the edge-weights can be integers, reals
45 /// or of other comparable numeric type. This algorithm is intended
46 /// to be used only for small values of \c k, since it is only
47 /// polynomial in k, not in the length of k (which is log k): in
48 /// order to find the minimum cost flow of value \c k it finds the
49 /// minimum cost flow of value \c i for every \c i between 0 and \c
52 ///\param Graph The directed graph type the algorithm runs on.
53 ///\param LengthMap The type of the length map.
54 ///\param CapacityMap The capacity map type.
56 ///\author Attila Bernath
57 template <typename Graph, typename LengthMap, typename CapacityMap>
60 typedef typename LengthMap::Value Length;
62 //Warning: this should be integer type
63 typedef typename CapacityMap::Value Capacity;
65 typedef typename Graph::Node Node;
66 typedef typename Graph::NodeIt NodeIt;
67 typedef typename Graph::Edge Edge;
68 typedef typename Graph::OutEdgeIt OutEdgeIt;
69 typedef typename Graph::template EdgeMap<int> EdgeIntMap;
71 typedef ResGraphAdaptor<const Graph,int,CapacityMap,EdgeIntMap> ResGW;
72 typedef typename ResGW::Edge ResGraphEdge;
77 const LengthMap& length;
78 const CapacityMap& capacity;
81 typedef typename Graph::template NodeMap<Length> PotentialMap;
82 PotentialMap potential;
90 typedef typename Graph::template NodeMap<Length> NodeMap;
92 const LengthMap &length;
95 typedef typename LengthMap::Key Key;
96 typedef typename LengthMap::Value Value;
98 ModLengthMap(const ResGW& _g,
99 const LengthMap &_length, const NodeMap &_pot) :
100 g(_g), /*rev(_rev),*/ length(_length), pot(_pot) { }
102 Value operator[](typename ResGW::Edge e) const {
104 return length[e]-(pot[g.target(e)]-pot[g.source(e)]);
106 return -length[e]-(pot[g.target(e)]-pot[g.source(e)]);
112 ModLengthMap mod_length;
113 Dijkstra<ResGW, ModLengthMap> dijkstra;
117 /*! \brief The constructor of the class.
119 \param _g The directed graph the algorithm runs on.
120 \param _length The length (weight or cost) of the edges.
121 \param _cap The capacity of the edges.
122 \param _s Source node.
123 \param _t Target node.
125 MinCostFlow(Graph& _g, LengthMap& _length, CapacityMap& _cap,
127 g(_g), length(_length), capacity(_cap), flow(_g), potential(_g),
129 res_graph(g, capacity, flow),
130 mod_length(res_graph, length, potential),
131 dijkstra(res_graph, mod_length) {
135 /*! Tries to augment the flow between s and t by 1.
136 The return value shows if the augmentation is successful.
140 if (!dijkstra.reached(t)) {
142 //Unsuccessful augmentation.
146 //We have to change the potential
147 for(typename ResGW::NodeIt n(res_graph); n!=INVALID; ++n)
148 potential.set(n, potential[n]+dijkstra.distMap()[n]);
150 //Augmenting on the shortest path
154 e = dijkstra.pred(n);
155 n = dijkstra.predNode(n);
156 res_graph.augment(e,1);
157 //Let's update the total length
158 if (res_graph.forward(e))
159 total_length += length[e];
161 total_length -= length[e];
168 /*! \brief Runs the algorithm.
171 Returns k if there is a flow of value at least k from s to t.
172 Otherwise it returns the maximum value of a flow from s to t.
174 \param k The value of the flow we are looking for.
176 \todo May be it does make sense to be able to start with a nonzero
177 feasible primal-dual solution pair as well.
179 \todo If the actual flow value is bigger than k, then everything is
180 cleared and the algorithm starts from zero flow. Is it a good approach?
183 if (flowValue()>k) reset();
184 while (flowValue()<k && augment()) { }
188 /*! \brief The class is reset to zero flow and potential.
189 The class is reset to zero flow and potential.
193 for (typename Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0);
194 for (typename Graph::NodeIt n(g); n!=INVALID; ++n) potential.set(n, 0);
197 /*! Returns the value of the actual flow.
199 int flowValue() const {
201 for (typename Graph::OutEdgeIt e(g, s); e!=INVALID; ++e) i+=flow[e];
202 for (typename Graph::InEdgeIt e(g, s); e!=INVALID; ++e) i-=flow[e];
206 /// Total weight of the found flow.
208 /// This function gives back the total weight of the found flow.
209 Length totalLength(){
213 ///Returns a const reference to the EdgeMap \c flow.
215 ///Returns a const reference to the EdgeMap \c flow.
216 const EdgeIntMap &getFlow() const { return flow;}
218 /*! \brief Returns a const reference to the NodeMap \c potential (the dual solution).
220 Returns a const reference to the NodeMap \c potential (the dual solution).
222 const PotentialMap &getPotential() const { return potential;}
224 /*! \brief Checking the complementary slackness optimality criteria.
226 This function checks, whether the given flow and potential
227 satisfy the complementary slackness conditions (i.e. these are optimal).
228 This function only checks optimality, doesn't bother with feasibility.
231 bool checkComplementarySlackness(){
234 for(typename Graph::EdgeIt e(g); e!=INVALID; ++e) {
236 mod_pot = length[e]-potential[g.target(e)]+potential[g.source(e)];
238 if (0<fl_e && fl_e<capacity[e]) {
239 /// \todo better comparison is needed for real types, moreover,
240 /// this comparison here is superfluous.
245 if (mod_pot > 0 && fl_e != 0)
247 if (mod_pot < 0 && fl_e != capacity[e])
254 }; //class MinCostFlow
260 #endif //LEMON_MIN_COST_FLOW_H