1.1 --- a/lemon/cycle_canceling.h Wed Feb 27 11:39:03 2008 +0000
1.2 +++ b/lemon/cycle_canceling.h Thu Feb 28 02:54:27 2008 +0000
1.3 @@ -52,11 +52,8 @@
1.4 /// \warning
1.5 /// - Edge capacities and costs should be \e non-negative \e integers.
1.6 /// - Supply values should be \e signed \e integers.
1.7 - /// - \c LowerMap::Value must be convertible to \c CapacityMap::Value.
1.8 - /// - \c CapacityMap::Value and \c SupplyMap::Value must be
1.9 - /// convertible to each other.
1.10 - /// - All value types must be convertible to \c CostMap::Value, which
1.11 - /// must be signed type.
1.12 + /// - The value types of the maps should be convertible to each other.
1.13 + /// - \c CostMap::Value must be signed type.
1.14 ///
1.15 /// \note By default the \ref BellmanFord "Bellman-Ford" algorithm is
1.16 /// used for negative cycle detection with limited iteration number.
1.17 @@ -94,6 +91,8 @@
1.18
1.19 /// The type of the flow map.
1.20 typedef typename Graph::template EdgeMap<Capacity> FlowMap;
1.21 + /// The type of the potential map.
1.22 + typedef typename Graph::template NodeMap<Cost> PotentialMap;
1.23
1.24 private:
1.25
1.26 @@ -143,18 +142,22 @@
1.27 bool _valid_supply;
1.28
1.29 // Edge map of the current flow
1.30 - FlowMap _flow;
1.31 + FlowMap *_flow;
1.32 + bool _local_flow;
1.33 + // Node map of the current potentials
1.34 + PotentialMap *_potential;
1.35 + bool _local_potential;
1.36
1.37 // The residual graph
1.38 - ResGraph _res_graph;
1.39 + ResGraph *_res_graph;
1.40 // The residual cost map
1.41 ResidualCostMap _res_cost;
1.42
1.43 public:
1.44
1.45 - /// \brief General constructor of the class (with lower bounds).
1.46 + /// \brief General constructor (with lower bounds).
1.47 ///
1.48 - /// General constructor of the class (with lower bounds).
1.49 + /// General constructor (with lower bounds).
1.50 ///
1.51 /// \param graph The directed graph the algorithm runs on.
1.52 /// \param lower The lower bounds of the edges.
1.53 @@ -167,8 +170,8 @@
1.54 const CostMap &cost,
1.55 const SupplyMap &supply ) :
1.56 _graph(graph), _lower(&lower), _capacity(graph), _cost(cost),
1.57 - _supply(graph), _flow(graph, 0),
1.58 - _res_graph(graph, _capacity, _flow), _res_cost(_cost)
1.59 + _supply(graph), _flow(0), _local_flow(false),
1.60 + _potential(0), _local_potential(false), _res_cost(_cost)
1.61 {
1.62 // Removing non-zero lower bounds
1.63 _capacity = subMap(capacity, lower);
1.64 @@ -184,9 +187,9 @@
1.65 _valid_supply = sum == 0;
1.66 }
1.67
1.68 - /// \brief General constructor of the class (without lower bounds).
1.69 + /// \brief General constructor (without lower bounds).
1.70 ///
1.71 - /// General constructor of the class (without lower bounds).
1.72 + /// General constructor (without lower bounds).
1.73 ///
1.74 /// \param graph The directed graph the algorithm runs on.
1.75 /// \param capacity The capacities (upper bounds) of the edges.
1.76 @@ -197,8 +200,8 @@
1.77 const CostMap &cost,
1.78 const SupplyMap &supply ) :
1.79 _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost),
1.80 - _supply(supply), _flow(graph, 0),
1.81 - _res_graph(graph, _capacity, _flow), _res_cost(_cost)
1.82 + _supply(supply), _flow(0), _local_flow(false),
1.83 + _potential(0), _local_potential(false), _res_cost(_cost)
1.84 {
1.85 // Checking the sum of supply values
1.86 Supply sum = 0;
1.87 @@ -206,9 +209,9 @@
1.88 _valid_supply = sum == 0;
1.89 }
1.90
1.91 - /// \brief Simple constructor of the class (with lower bounds).
1.92 + /// \brief Simple constructor (with lower bounds).
1.93 ///
1.94 - /// Simple constructor of the class (with lower bounds).
1.95 + /// Simple constructor (with lower bounds).
1.96 ///
1.97 /// \param graph The directed graph the algorithm runs on.
1.98 /// \param lower The lower bounds of the edges.
1.99 @@ -225,8 +228,8 @@
1.100 Node s, Node t,
1.101 Supply flow_value ) :
1.102 _graph(graph), _lower(&lower), _capacity(graph), _cost(cost),
1.103 - _supply(graph), _flow(graph, 0),
1.104 - _res_graph(graph, _capacity, _flow), _res_cost(_cost)
1.105 + _supply(graph), _flow(0), _local_flow(false),
1.106 + _potential(0), _local_potential(false), _res_cost(_cost)
1.107 {
1.108 // Removing non-zero lower bounds
1.109 _capacity = subMap(capacity, lower);
1.110 @@ -243,9 +246,9 @@
1.111 _valid_supply = true;
1.112 }
1.113
1.114 - /// \brief Simple constructor of the class (without lower bounds).
1.115 + /// \brief Simple constructor (without lower bounds).
1.116 ///
1.117 - /// Simple constructor of the class (without lower bounds).
1.118 + /// Simple constructor (without lower bounds).
1.119 ///
1.120 /// \param graph The directed graph the algorithm runs on.
1.121 /// \param capacity The capacities (upper bounds) of the edges.
1.122 @@ -260,14 +263,55 @@
1.123 Node s, Node t,
1.124 Supply flow_value ) :
1.125 _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost),
1.126 - _supply(graph, 0), _flow(graph, 0),
1.127 - _res_graph(graph, _capacity, _flow), _res_cost(_cost)
1.128 + _supply(graph, 0), _flow(0), _local_flow(false),
1.129 + _potential(0), _local_potential(false), _res_cost(_cost)
1.130 {
1.131 _supply[s] = flow_value;
1.132 _supply[t] = -flow_value;
1.133 _valid_supply = true;
1.134 }
1.135
1.136 + /// Destructor.
1.137 + ~CycleCanceling() {
1.138 + if (_local_flow) delete _flow;
1.139 + if (_local_potential) delete _potential;
1.140 + delete _res_graph;
1.141 + }
1.142 +
1.143 + /// \brief Sets the flow map.
1.144 + ///
1.145 + /// Sets the flow map.
1.146 + ///
1.147 + /// \return \c (*this)
1.148 + CycleCanceling& flowMap(FlowMap &map) {
1.149 + if (_local_flow) {
1.150 + delete _flow;
1.151 + _local_flow = false;
1.152 + }
1.153 + _flow = ↦
1.154 + return *this;
1.155 + }
1.156 +
1.157 + /// \brief Sets the potential map.
1.158 + ///
1.159 + /// Sets the potential map.
1.160 + ///
1.161 + /// \return \c (*this)
1.162 + CycleCanceling& potentialMap(PotentialMap &map) {
1.163 + if (_local_potential) {
1.164 + delete _potential;
1.165 + _local_potential = false;
1.166 + }
1.167 + _potential = ↦
1.168 + return *this;
1.169 + }
1.170 +
1.171 + /// \name Execution control
1.172 + /// The only way to execute the algorithm is to call the run()
1.173 + /// function.
1.174 +
1.175 + /// @{
1.176 +
1.177 /// \brief Runs the algorithm.
1.178 ///
1.179 /// Runs the algorithm.
1.180 @@ -281,6 +325,15 @@
1.181 return init() && start(min_mean_cc);
1.182 }
1.183
1.184 + /// @}
1.185 +
1.186 + /// \name Query Functions
1.187 + /// The result of the algorithm can be obtained using these
1.188 + /// functions.
1.189 + /// \n run() must be called before using them.
1.190 +
1.191 + /// @{
1.192 +
1.193 /// \brief Returns a const reference to the edge map storing the
1.194 /// found flow.
1.195 ///
1.196 @@ -288,7 +341,36 @@
1.197 ///
1.198 /// \pre \ref run() must be called before using this function.
1.199 const FlowMap& flowMap() const {
1.200 - return _flow;
1.201 + return *_flow;
1.202 + }
1.203 +
1.204 + /// \brief Returns a const reference to the node map storing the
1.205 + /// found potentials (the dual solution).
1.206 + ///
1.207 + /// Returns a const reference to the node map storing the found
1.208 + /// potentials (the dual solution).
1.209 + ///
1.210 + /// \pre \ref run() must be called before using this function.
1.211 + const PotentialMap& potentialMap() const {
1.212 + return *_potential;
1.213 + }
1.214 +
1.215 + /// \brief Returns the flow on the edge.
1.216 + ///
1.217 + /// Returns the flow on the edge.
1.218 + ///
1.219 + /// \pre \ref run() must be called before using this function.
1.220 + Capacity flow(const Edge& edge) const {
1.221 + return (*_flow)[edge];
1.222 + }
1.223 +
1.224 + /// \brief Returns the potential of the node.
1.225 + ///
1.226 + /// Returns the potential of the node.
1.227 + ///
1.228 + /// \pre \ref run() must be called before using this function.
1.229 + Cost potential(const Node& node) const {
1.230 + return (*_potential)[node];
1.231 }
1.232
1.233 /// \brief Returns the total cost of the found flow.
1.234 @@ -300,29 +382,50 @@
1.235 Cost totalCost() const {
1.236 Cost c = 0;
1.237 for (EdgeIt e(_graph); e != INVALID; ++e)
1.238 - c += _flow[e] * _cost[e];
1.239 + c += (*_flow)[e] * _cost[e];
1.240 return c;
1.241 }
1.242
1.243 + /// @}
1.244 +
1.245 private:
1.246
1.247 /// Initializes the algorithm.
1.248 bool init() {
1.249 if (!_valid_supply) return false;
1.250
1.251 + // Initializing flow and potential maps
1.252 + if (!_flow) {
1.253 + _flow = new FlowMap(_graph);
1.254 + _local_flow = true;
1.255 + }
1.256 + if (!_potential) {
1.257 + _potential = new PotentialMap(_graph);
1.258 + _local_potential = true;
1.259 + }
1.260 +
1.261 + _res_graph = new ResGraph(_graph, _capacity, *_flow);
1.262 +
1.263 // Finding a feasible flow using Circulation
1.264 Circulation< Graph, ConstMap<Edge, Capacity>, CapacityEdgeMap,
1.265 SupplyMap >
1.266 - circulation( _graph, constMap<Edge>((Capacity)0), _capacity,
1.267 + circulation( _graph, constMap<Edge>(Capacity(0)), _capacity,
1.268 _supply );
1.269 - return circulation.flowMap(_flow).run();
1.270 + return circulation.flowMap(*_flow).run();
1.271 }
1.272
1.273 bool start(bool min_mean_cc) {
1.274 if (min_mean_cc)
1.275 - return startMinMean();
1.276 + startMinMean();
1.277 else
1.278 - return start();
1.279 + start();
1.280 +
1.281 + // Handling non-zero lower bounds
1.282 + if (_lower) {
1.283 + for (EdgeIt e(_graph); e != INVALID; ++e)
1.284 + (*_flow)[e] += (*_lower)[e];
1.285 + }
1.286 + return true;
1.287 }
1.288
1.289 /// \brief Executes the algorithm using \ref BellmanFord.
1.290 @@ -330,16 +433,16 @@
1.291 /// Executes the algorithm using the \ref BellmanFord
1.292 /// "Bellman-Ford" algorithm for negative cycle detection with
1.293 /// successively larger limit for the number of iterations.
1.294 - bool start() {
1.295 - typename BellmanFord<ResGraph, ResidualCostMap>::PredMap pred(_res_graph);
1.296 - typename ResGraph::template NodeMap<int> visited(_res_graph);
1.297 + void start() {
1.298 + typename BellmanFord<ResGraph, ResidualCostMap>::PredMap pred(*_res_graph);
1.299 + typename ResGraph::template NodeMap<int> visited(*_res_graph);
1.300 std::vector<ResEdge> cycle;
1.301 int node_num = countNodes(_graph);
1.302
1.303 int length_bound = BF_FIRST_LIMIT;
1.304 bool optimal = false;
1.305 while (!optimal) {
1.306 - BellmanFord<ResGraph, ResidualCostMap> bf(_res_graph, _res_cost);
1.307 + BellmanFord<ResGraph, ResidualCostMap> bf(*_res_graph, _res_cost);
1.308 bf.predMap(pred);
1.309 bf.init(0);
1.310 int iter_num = 0;
1.311 @@ -356,22 +459,26 @@
1.312 }
1.313 }
1.314 if (real_iter_num < curr_iter_num) {
1.315 + // Optimal flow is found
1.316 optimal = true;
1.317 + // Setting node potentials
1.318 + for (NodeIt n(_graph); n != INVALID; ++n)
1.319 + (*_potential)[n] = bf.dist(n);
1.320 break;
1.321 } else {
1.322 // Searching for node disjoint negative cycles
1.323 - for (ResNodeIt n(_res_graph); n != INVALID; ++n)
1.324 + for (ResNodeIt n(*_res_graph); n != INVALID; ++n)
1.325 visited[n] = 0;
1.326 int id = 0;
1.327 - for (ResNodeIt n(_res_graph); n != INVALID; ++n) {
1.328 + for (ResNodeIt n(*_res_graph); n != INVALID; ++n) {
1.329 if (visited[n] > 0) continue;
1.330 visited[n] = ++id;
1.331 ResNode u = pred[n] == INVALID ?
1.332 - INVALID : _res_graph.source(pred[n]);
1.333 + INVALID : _res_graph->source(pred[n]);
1.334 while (u != INVALID && visited[u] == 0) {
1.335 visited[u] = id;
1.336 u = pred[u] == INVALID ?
1.337 - INVALID : _res_graph.source(pred[u]);
1.338 + INVALID : _res_graph->source(pred[u]);
1.339 }
1.340 if (u != INVALID && visited[u] == id) {
1.341 // Finding the negative cycle
1.342 @@ -379,16 +486,16 @@
1.343 cycle.clear();
1.344 ResEdge e = pred[u];
1.345 cycle.push_back(e);
1.346 - Capacity d = _res_graph.rescap(e);
1.347 - while (_res_graph.source(e) != u) {
1.348 - cycle.push_back(e = pred[_res_graph.source(e)]);
1.349 - if (_res_graph.rescap(e) < d)
1.350 - d = _res_graph.rescap(e);
1.351 + Capacity d = _res_graph->rescap(e);
1.352 + while (_res_graph->source(e) != u) {
1.353 + cycle.push_back(e = pred[_res_graph->source(e)]);
1.354 + if (_res_graph->rescap(e) < d)
1.355 + d = _res_graph->rescap(e);
1.356 }
1.357
1.358 // Augmenting along the cycle
1.359 for (int i = 0; i < int(cycle.size()); ++i)
1.360 - _res_graph.augment(cycle[i], d);
1.361 + _res_graph->augment(cycle[i], d);
1.362 }
1.363 }
1.364 }
1.365 @@ -397,22 +504,15 @@
1.366 length_bound = int(length_bound * BF_ALPHA);
1.367 }
1.368 }
1.369 -
1.370 - // Handling non-zero lower bounds
1.371 - if (_lower) {
1.372 - for (EdgeIt e(_graph); e != INVALID; ++e)
1.373 - _flow[e] += (*_lower)[e];
1.374 - }
1.375 - return true;
1.376 }
1.377
1.378 /// \brief Executes the algorithm using \ref MinMeanCycle.
1.379 ///
1.380 /// Executes the algorithm using \ref MinMeanCycle for negative
1.381 /// cycle detection.
1.382 - bool startMinMean() {
1.383 + void startMinMean() {
1.384 typedef Path<ResGraph> ResPath;
1.385 - MinMeanCycle<ResGraph, ResidualCostMap> mmc(_res_graph, _res_cost);
1.386 + MinMeanCycle<ResGraph, ResidualCostMap> mmc(*_res_graph, _res_cost);
1.387 ResPath cycle;
1.388
1.389 mmc.cyclePath(cycle).init();
1.390 @@ -425,13 +525,13 @@
1.391 // along the cycle
1.392 Capacity delta = 0;
1.393 for (typename ResPath::EdgeIt e(cycle); e != INVALID; ++e) {
1.394 - if (delta == 0 || _res_graph.rescap(e) < delta)
1.395 - delta = _res_graph.rescap(e);
1.396 + if (delta == 0 || _res_graph->rescap(e) < delta)
1.397 + delta = _res_graph->rescap(e);
1.398 }
1.399
1.400 // Augmenting along the cycle
1.401 for (typename ResPath::EdgeIt e(cycle); e != INVALID; ++e)
1.402 - _res_graph.augment(e, delta);
1.403 + _res_graph->augment(e, delta);
1.404
1.405 // Finding the minimum cycle mean for the modified residual
1.406 // graph
1.407 @@ -440,12 +540,11 @@
1.408 }
1.409 }
1.410
1.411 - // Handling non-zero lower bounds
1.412 - if (_lower) {
1.413 - for (EdgeIt e(_graph); e != INVALID; ++e)
1.414 - _flow[e] += (*_lower)[e];
1.415 - }
1.416 - return true;
1.417 + // Computing node potentials
1.418 + BellmanFord<ResGraph, ResidualCostMap> bf(*_res_graph, _res_cost);
1.419 + bf.init(0); bf.start();
1.420 + for (NodeIt n(_graph); n != INVALID; ++n)
1.421 + (*_potential)[n] = bf.dist(n);
1.422 }
1.423
1.424 }; //class CycleCanceling