kpeter@814
|
1 |
/* -*- C++ -*-
|
kpeter@814
|
2 |
*
|
kpeter@814
|
3 |
* This file is a part of LEMON, a generic C++ optimization library
|
kpeter@814
|
4 |
*
|
kpeter@814
|
5 |
* Copyright (C) 2003-2008
|
kpeter@814
|
6 |
* Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
|
kpeter@814
|
7 |
* (Egervary Research Group on Combinatorial Optimization, EGRES).
|
kpeter@814
|
8 |
*
|
kpeter@814
|
9 |
* Permission to use, modify and distribute this software is granted
|
kpeter@814
|
10 |
* provided that this copyright notice appears in all copies. For
|
kpeter@814
|
11 |
* precise terms see the accompanying LICENSE file.
|
kpeter@814
|
12 |
*
|
kpeter@814
|
13 |
* This software is provided "AS IS" with no warranty of any kind,
|
kpeter@814
|
14 |
* express or implied, and with no claim as to its suitability for any
|
kpeter@814
|
15 |
* purpose.
|
kpeter@814
|
16 |
*
|
kpeter@814
|
17 |
*/
|
kpeter@814
|
18 |
|
kpeter@814
|
19 |
#ifndef LEMON_CYCLE_CANCELING_H
|
kpeter@814
|
20 |
#define LEMON_CYCLE_CANCELING_H
|
kpeter@814
|
21 |
|
kpeter@814
|
22 |
/// \ingroup min_cost_flow
|
kpeter@814
|
23 |
///
|
kpeter@814
|
24 |
/// \file
|
kpeter@814
|
25 |
/// \brief Cycle-canceling algorithm for finding a minimum cost flow.
|
kpeter@814
|
26 |
|
kpeter@814
|
27 |
#include <vector>
|
kpeter@814
|
28 |
#include <lemon/adaptors.h>
|
kpeter@814
|
29 |
#include <lemon/path.h>
|
kpeter@814
|
30 |
|
kpeter@814
|
31 |
#include <lemon/circulation.h>
|
kpeter@814
|
32 |
#include <lemon/bellman_ford.h>
|
kpeter@814
|
33 |
#include <lemon/howard.h>
|
kpeter@814
|
34 |
|
kpeter@814
|
35 |
namespace lemon {
|
kpeter@814
|
36 |
|
kpeter@814
|
37 |
/// \addtogroup min_cost_flow
|
kpeter@814
|
38 |
/// @{
|
kpeter@814
|
39 |
|
kpeter@814
|
40 |
/// \brief Implementation of a cycle-canceling algorithm for
|
kpeter@814
|
41 |
/// finding a minimum cost flow.
|
kpeter@814
|
42 |
///
|
kpeter@814
|
43 |
/// \ref CycleCanceling implements a cycle-canceling algorithm for
|
kpeter@814
|
44 |
/// finding a minimum cost flow.
|
kpeter@814
|
45 |
///
|
kpeter@814
|
46 |
/// \tparam Digraph The digraph type the algorithm runs on.
|
kpeter@814
|
47 |
/// \tparam LowerMap The type of the lower bound map.
|
kpeter@814
|
48 |
/// \tparam CapacityMap The type of the capacity (upper bound) map.
|
kpeter@814
|
49 |
/// \tparam CostMap The type of the cost (length) map.
|
kpeter@814
|
50 |
/// \tparam SupplyMap The type of the supply map.
|
kpeter@814
|
51 |
///
|
kpeter@814
|
52 |
/// \warning
|
kpeter@814
|
53 |
/// - Arc capacities and costs should be \e non-negative \e integers.
|
kpeter@814
|
54 |
/// - Supply values should be \e signed \e integers.
|
kpeter@814
|
55 |
/// - The value types of the maps should be convertible to each other.
|
kpeter@814
|
56 |
/// - \c CostMap::Value must be signed type.
|
kpeter@814
|
57 |
///
|
kpeter@814
|
58 |
/// \note By default the \ref BellmanFord "Bellman-Ford" algorithm is
|
kpeter@814
|
59 |
/// used for negative cycle detection with limited iteration number.
|
kpeter@814
|
60 |
/// However \ref CycleCanceling also provides the "Minimum Mean
|
kpeter@814
|
61 |
/// Cycle-Canceling" algorithm, which is \e strongly \e polynomial,
|
kpeter@814
|
62 |
/// but rather slower in practice.
|
kpeter@814
|
63 |
/// To use this version of the algorithm, call \ref run() with \c true
|
kpeter@814
|
64 |
/// parameter.
|
kpeter@814
|
65 |
///
|
kpeter@814
|
66 |
/// \author Peter Kovacs
|
kpeter@814
|
67 |
template < typename Digraph,
|
kpeter@814
|
68 |
typename LowerMap = typename Digraph::template ArcMap<int>,
|
kpeter@814
|
69 |
typename CapacityMap = typename Digraph::template ArcMap<int>,
|
kpeter@814
|
70 |
typename CostMap = typename Digraph::template ArcMap<int>,
|
kpeter@814
|
71 |
typename SupplyMap = typename Digraph::template NodeMap<int> >
|
kpeter@814
|
72 |
class CycleCanceling
|
kpeter@814
|
73 |
{
|
kpeter@814
|
74 |
TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
|
kpeter@814
|
75 |
|
kpeter@814
|
76 |
typedef typename CapacityMap::Value Capacity;
|
kpeter@814
|
77 |
typedef typename CostMap::Value Cost;
|
kpeter@814
|
78 |
typedef typename SupplyMap::Value Supply;
|
kpeter@814
|
79 |
typedef typename Digraph::template ArcMap<Capacity> CapacityArcMap;
|
kpeter@814
|
80 |
typedef typename Digraph::template NodeMap<Supply> SupplyNodeMap;
|
kpeter@814
|
81 |
|
kpeter@814
|
82 |
typedef ResidualDigraph< const Digraph,
|
kpeter@814
|
83 |
CapacityArcMap, CapacityArcMap > ResDigraph;
|
kpeter@814
|
84 |
typedef typename ResDigraph::Node ResNode;
|
kpeter@814
|
85 |
typedef typename ResDigraph::NodeIt ResNodeIt;
|
kpeter@814
|
86 |
typedef typename ResDigraph::Arc ResArc;
|
kpeter@814
|
87 |
typedef typename ResDigraph::ArcIt ResArcIt;
|
kpeter@814
|
88 |
|
kpeter@814
|
89 |
public:
|
kpeter@814
|
90 |
|
kpeter@814
|
91 |
/// The type of the flow map.
|
kpeter@814
|
92 |
typedef typename Digraph::template ArcMap<Capacity> FlowMap;
|
kpeter@814
|
93 |
/// The type of the potential map.
|
kpeter@814
|
94 |
typedef typename Digraph::template NodeMap<Cost> PotentialMap;
|
kpeter@814
|
95 |
|
kpeter@814
|
96 |
private:
|
kpeter@814
|
97 |
|
kpeter@814
|
98 |
/// \brief Map adaptor class for handling residual arc costs.
|
kpeter@814
|
99 |
///
|
kpeter@814
|
100 |
/// Map adaptor class for handling residual arc costs.
|
kpeter@814
|
101 |
class ResidualCostMap : public MapBase<ResArc, Cost>
|
kpeter@814
|
102 |
{
|
kpeter@814
|
103 |
private:
|
kpeter@814
|
104 |
|
kpeter@814
|
105 |
const CostMap &_cost_map;
|
kpeter@814
|
106 |
|
kpeter@814
|
107 |
public:
|
kpeter@814
|
108 |
|
kpeter@814
|
109 |
///\e
|
kpeter@814
|
110 |
ResidualCostMap(const CostMap &cost_map) : _cost_map(cost_map) {}
|
kpeter@814
|
111 |
|
kpeter@814
|
112 |
///\e
|
kpeter@814
|
113 |
Cost operator[](const ResArc &e) const {
|
kpeter@814
|
114 |
return ResDigraph::forward(e) ? _cost_map[e] : -_cost_map[e];
|
kpeter@814
|
115 |
}
|
kpeter@814
|
116 |
|
kpeter@814
|
117 |
}; //class ResidualCostMap
|
kpeter@814
|
118 |
|
kpeter@814
|
119 |
private:
|
kpeter@814
|
120 |
|
kpeter@814
|
121 |
// The maximum number of iterations for the first execution of the
|
kpeter@814
|
122 |
// Bellman-Ford algorithm. It should be at least 2.
|
kpeter@814
|
123 |
static const int BF_FIRST_LIMIT = 2;
|
kpeter@814
|
124 |
// The iteration limit for the Bellman-Ford algorithm is multiplied
|
kpeter@814
|
125 |
// by BF_LIMIT_FACTOR/100 in every round.
|
kpeter@814
|
126 |
static const int BF_LIMIT_FACTOR = 150;
|
kpeter@814
|
127 |
|
kpeter@814
|
128 |
private:
|
kpeter@814
|
129 |
|
kpeter@814
|
130 |
// The digraph the algorithm runs on
|
kpeter@814
|
131 |
const Digraph &_graph;
|
kpeter@814
|
132 |
// The original lower bound map
|
kpeter@814
|
133 |
const LowerMap *_lower;
|
kpeter@814
|
134 |
// The modified capacity map
|
kpeter@814
|
135 |
CapacityArcMap _capacity;
|
kpeter@814
|
136 |
// The original cost map
|
kpeter@814
|
137 |
const CostMap &_cost;
|
kpeter@814
|
138 |
// The modified supply map
|
kpeter@814
|
139 |
SupplyNodeMap _supply;
|
kpeter@814
|
140 |
bool _valid_supply;
|
kpeter@814
|
141 |
|
kpeter@814
|
142 |
// Arc map of the current flow
|
kpeter@814
|
143 |
FlowMap *_flow;
|
kpeter@814
|
144 |
bool _local_flow;
|
kpeter@814
|
145 |
// Node map of the current potentials
|
kpeter@814
|
146 |
PotentialMap *_potential;
|
kpeter@814
|
147 |
bool _local_potential;
|
kpeter@814
|
148 |
|
kpeter@814
|
149 |
// The residual digraph
|
kpeter@814
|
150 |
ResDigraph *_res_graph;
|
kpeter@814
|
151 |
// The residual cost map
|
kpeter@814
|
152 |
ResidualCostMap _res_cost;
|
kpeter@814
|
153 |
|
kpeter@814
|
154 |
public:
|
kpeter@814
|
155 |
|
kpeter@814
|
156 |
/// \brief General constructor (with lower bounds).
|
kpeter@814
|
157 |
///
|
kpeter@814
|
158 |
/// General constructor (with lower bounds).
|
kpeter@814
|
159 |
///
|
kpeter@814
|
160 |
/// \param digraph The digraph the algorithm runs on.
|
kpeter@814
|
161 |
/// \param lower The lower bounds of the arcs.
|
kpeter@814
|
162 |
/// \param capacity The capacities (upper bounds) of the arcs.
|
kpeter@814
|
163 |
/// \param cost The cost (length) values of the arcs.
|
kpeter@814
|
164 |
/// \param supply The supply values of the nodes (signed).
|
kpeter@814
|
165 |
CycleCanceling( const Digraph &digraph,
|
kpeter@814
|
166 |
const LowerMap &lower,
|
kpeter@814
|
167 |
const CapacityMap &capacity,
|
kpeter@814
|
168 |
const CostMap &cost,
|
kpeter@814
|
169 |
const SupplyMap &supply ) :
|
kpeter@814
|
170 |
_graph(digraph), _lower(&lower), _capacity(digraph), _cost(cost),
|
kpeter@814
|
171 |
_supply(digraph), _flow(NULL), _local_flow(false),
|
kpeter@814
|
172 |
_potential(NULL), _local_potential(false),
|
kpeter@814
|
173 |
_res_graph(NULL), _res_cost(_cost)
|
kpeter@814
|
174 |
{
|
kpeter@814
|
175 |
// Check the sum of supply values
|
kpeter@814
|
176 |
Supply sum = 0;
|
kpeter@814
|
177 |
for (NodeIt n(_graph); n != INVALID; ++n) {
|
kpeter@814
|
178 |
_supply[n] = supply[n];
|
kpeter@814
|
179 |
sum += _supply[n];
|
kpeter@814
|
180 |
}
|
kpeter@814
|
181 |
_valid_supply = sum == 0;
|
kpeter@814
|
182 |
|
kpeter@814
|
183 |
// Remove non-zero lower bounds
|
kpeter@814
|
184 |
for (ArcIt e(_graph); e != INVALID; ++e) {
|
kpeter@814
|
185 |
_capacity[e] = capacity[e];
|
kpeter@814
|
186 |
if (lower[e] != 0) {
|
kpeter@814
|
187 |
_capacity[e] -= lower[e];
|
kpeter@814
|
188 |
_supply[_graph.source(e)] -= lower[e];
|
kpeter@814
|
189 |
_supply[_graph.target(e)] += lower[e];
|
kpeter@814
|
190 |
}
|
kpeter@814
|
191 |
}
|
kpeter@814
|
192 |
}
|
kpeter@814
|
193 |
/*
|
kpeter@814
|
194 |
/// \brief General constructor (without lower bounds).
|
kpeter@814
|
195 |
///
|
kpeter@814
|
196 |
/// General constructor (without lower bounds).
|
kpeter@814
|
197 |
///
|
kpeter@814
|
198 |
/// \param digraph The digraph the algorithm runs on.
|
kpeter@814
|
199 |
/// \param capacity The capacities (upper bounds) of the arcs.
|
kpeter@814
|
200 |
/// \param cost The cost (length) values of the arcs.
|
kpeter@814
|
201 |
/// \param supply The supply values of the nodes (signed).
|
kpeter@814
|
202 |
CycleCanceling( const Digraph &digraph,
|
kpeter@814
|
203 |
const CapacityMap &capacity,
|
kpeter@814
|
204 |
const CostMap &cost,
|
kpeter@814
|
205 |
const SupplyMap &supply ) :
|
kpeter@814
|
206 |
_graph(digraph), _lower(NULL), _capacity(capacity), _cost(cost),
|
kpeter@814
|
207 |
_supply(supply), _flow(NULL), _local_flow(false),
|
kpeter@814
|
208 |
_potential(NULL), _local_potential(false), _res_graph(NULL),
|
kpeter@814
|
209 |
_res_cost(_cost)
|
kpeter@814
|
210 |
{
|
kpeter@814
|
211 |
// Check the sum of supply values
|
kpeter@814
|
212 |
Supply sum = 0;
|
kpeter@814
|
213 |
for (NodeIt n(_graph); n != INVALID; ++n) sum += _supply[n];
|
kpeter@814
|
214 |
_valid_supply = sum == 0;
|
kpeter@814
|
215 |
}
|
kpeter@814
|
216 |
|
kpeter@814
|
217 |
/// \brief Simple constructor (with lower bounds).
|
kpeter@814
|
218 |
///
|
kpeter@814
|
219 |
/// Simple constructor (with lower bounds).
|
kpeter@814
|
220 |
///
|
kpeter@814
|
221 |
/// \param digraph The digraph the algorithm runs on.
|
kpeter@814
|
222 |
/// \param lower The lower bounds of the arcs.
|
kpeter@814
|
223 |
/// \param capacity The capacities (upper bounds) of the arcs.
|
kpeter@814
|
224 |
/// \param cost The cost (length) values of the arcs.
|
kpeter@814
|
225 |
/// \param s The source node.
|
kpeter@814
|
226 |
/// \param t The target node.
|
kpeter@814
|
227 |
/// \param flow_value The required amount of flow from node \c s
|
kpeter@814
|
228 |
/// to node \c t (i.e. the supply of \c s and the demand of \c t).
|
kpeter@814
|
229 |
CycleCanceling( const Digraph &digraph,
|
kpeter@814
|
230 |
const LowerMap &lower,
|
kpeter@814
|
231 |
const CapacityMap &capacity,
|
kpeter@814
|
232 |
const CostMap &cost,
|
kpeter@814
|
233 |
Node s, Node t,
|
kpeter@814
|
234 |
Supply flow_value ) :
|
kpeter@814
|
235 |
_graph(digraph), _lower(&lower), _capacity(capacity), _cost(cost),
|
kpeter@814
|
236 |
_supply(digraph, 0), _flow(NULL), _local_flow(false),
|
kpeter@814
|
237 |
_potential(NULL), _local_potential(false), _res_graph(NULL),
|
kpeter@814
|
238 |
_res_cost(_cost)
|
kpeter@814
|
239 |
{
|
kpeter@814
|
240 |
// Remove non-zero lower bounds
|
kpeter@814
|
241 |
_supply[s] = flow_value;
|
kpeter@814
|
242 |
_supply[t] = -flow_value;
|
kpeter@814
|
243 |
for (ArcIt e(_graph); e != INVALID; ++e) {
|
kpeter@814
|
244 |
if (lower[e] != 0) {
|
kpeter@814
|
245 |
_capacity[e] -= lower[e];
|
kpeter@814
|
246 |
_supply[_graph.source(e)] -= lower[e];
|
kpeter@814
|
247 |
_supply[_graph.target(e)] += lower[e];
|
kpeter@814
|
248 |
}
|
kpeter@814
|
249 |
}
|
kpeter@814
|
250 |
_valid_supply = true;
|
kpeter@814
|
251 |
}
|
kpeter@814
|
252 |
|
kpeter@814
|
253 |
/// \brief Simple constructor (without lower bounds).
|
kpeter@814
|
254 |
///
|
kpeter@814
|
255 |
/// Simple constructor (without lower bounds).
|
kpeter@814
|
256 |
///
|
kpeter@814
|
257 |
/// \param digraph The digraph the algorithm runs on.
|
kpeter@814
|
258 |
/// \param capacity The capacities (upper bounds) of the arcs.
|
kpeter@814
|
259 |
/// \param cost The cost (length) values of the arcs.
|
kpeter@814
|
260 |
/// \param s The source node.
|
kpeter@814
|
261 |
/// \param t The target node.
|
kpeter@814
|
262 |
/// \param flow_value The required amount of flow from node \c s
|
kpeter@814
|
263 |
/// to node \c t (i.e. the supply of \c s and the demand of \c t).
|
kpeter@814
|
264 |
CycleCanceling( const Digraph &digraph,
|
kpeter@814
|
265 |
const CapacityMap &capacity,
|
kpeter@814
|
266 |
const CostMap &cost,
|
kpeter@814
|
267 |
Node s, Node t,
|
kpeter@814
|
268 |
Supply flow_value ) :
|
kpeter@814
|
269 |
_graph(digraph), _lower(NULL), _capacity(capacity), _cost(cost),
|
kpeter@814
|
270 |
_supply(digraph, 0), _flow(NULL), _local_flow(false),
|
kpeter@814
|
271 |
_potential(NULL), _local_potential(false), _res_graph(NULL),
|
kpeter@814
|
272 |
_res_cost(_cost)
|
kpeter@814
|
273 |
{
|
kpeter@814
|
274 |
_supply[s] = flow_value;
|
kpeter@814
|
275 |
_supply[t] = -flow_value;
|
kpeter@814
|
276 |
_valid_supply = true;
|
kpeter@814
|
277 |
}
|
kpeter@814
|
278 |
*/
|
kpeter@814
|
279 |
/// Destructor.
|
kpeter@814
|
280 |
~CycleCanceling() {
|
kpeter@814
|
281 |
if (_local_flow) delete _flow;
|
kpeter@814
|
282 |
if (_local_potential) delete _potential;
|
kpeter@814
|
283 |
delete _res_graph;
|
kpeter@814
|
284 |
}
|
kpeter@814
|
285 |
|
kpeter@814
|
286 |
/// \brief Set the flow map.
|
kpeter@814
|
287 |
///
|
kpeter@814
|
288 |
/// Set the flow map.
|
kpeter@814
|
289 |
///
|
kpeter@814
|
290 |
/// \return \c (*this)
|
kpeter@814
|
291 |
CycleCanceling& flowMap(FlowMap &map) {
|
kpeter@814
|
292 |
if (_local_flow) {
|
kpeter@814
|
293 |
delete _flow;
|
kpeter@814
|
294 |
_local_flow = false;
|
kpeter@814
|
295 |
}
|
kpeter@814
|
296 |
_flow = ↦
|
kpeter@814
|
297 |
return *this;
|
kpeter@814
|
298 |
}
|
kpeter@814
|
299 |
|
kpeter@814
|
300 |
/// \brief Set the potential map.
|
kpeter@814
|
301 |
///
|
kpeter@814
|
302 |
/// Set the potential map.
|
kpeter@814
|
303 |
///
|
kpeter@814
|
304 |
/// \return \c (*this)
|
kpeter@814
|
305 |
CycleCanceling& potentialMap(PotentialMap &map) {
|
kpeter@814
|
306 |
if (_local_potential) {
|
kpeter@814
|
307 |
delete _potential;
|
kpeter@814
|
308 |
_local_potential = false;
|
kpeter@814
|
309 |
}
|
kpeter@814
|
310 |
_potential = ↦
|
kpeter@814
|
311 |
return *this;
|
kpeter@814
|
312 |
}
|
kpeter@814
|
313 |
|
kpeter@814
|
314 |
/// \name Execution control
|
kpeter@814
|
315 |
|
kpeter@814
|
316 |
/// @{
|
kpeter@814
|
317 |
|
kpeter@814
|
318 |
/// \brief Run the algorithm.
|
kpeter@814
|
319 |
///
|
kpeter@814
|
320 |
/// Run the algorithm.
|
kpeter@814
|
321 |
///
|
kpeter@814
|
322 |
/// \param min_mean_cc Set this parameter to \c true to run the
|
kpeter@814
|
323 |
/// "Minimum Mean Cycle-Canceling" algorithm, which is strongly
|
kpeter@814
|
324 |
/// polynomial, but rather slower in practice.
|
kpeter@814
|
325 |
///
|
kpeter@814
|
326 |
/// \return \c true if a feasible flow can be found.
|
kpeter@814
|
327 |
bool run(bool min_mean_cc = false) {
|
kpeter@814
|
328 |
return init() && start(min_mean_cc);
|
kpeter@814
|
329 |
}
|
kpeter@814
|
330 |
|
kpeter@814
|
331 |
/// @}
|
kpeter@814
|
332 |
|
kpeter@814
|
333 |
/// \name Query Functions
|
kpeter@814
|
334 |
/// The result of the algorithm can be obtained using these
|
kpeter@814
|
335 |
/// functions.\n
|
kpeter@814
|
336 |
/// \ref lemon::CycleCanceling::run() "run()" must be called before
|
kpeter@814
|
337 |
/// using them.
|
kpeter@814
|
338 |
|
kpeter@814
|
339 |
/// @{
|
kpeter@814
|
340 |
|
kpeter@814
|
341 |
/// \brief Return a const reference to the arc map storing the
|
kpeter@814
|
342 |
/// found flow.
|
kpeter@814
|
343 |
///
|
kpeter@814
|
344 |
/// Return a const reference to the arc map storing the found flow.
|
kpeter@814
|
345 |
///
|
kpeter@814
|
346 |
/// \pre \ref run() must be called before using this function.
|
kpeter@814
|
347 |
const FlowMap& flowMap() const {
|
kpeter@814
|
348 |
return *_flow;
|
kpeter@814
|
349 |
}
|
kpeter@814
|
350 |
|
kpeter@814
|
351 |
/// \brief Return a const reference to the node map storing the
|
kpeter@814
|
352 |
/// found potentials (the dual solution).
|
kpeter@814
|
353 |
///
|
kpeter@814
|
354 |
/// Return a const reference to the node map storing the found
|
kpeter@814
|
355 |
/// potentials (the dual solution).
|
kpeter@814
|
356 |
///
|
kpeter@814
|
357 |
/// \pre \ref run() must be called before using this function.
|
kpeter@814
|
358 |
const PotentialMap& potentialMap() const {
|
kpeter@814
|
359 |
return *_potential;
|
kpeter@814
|
360 |
}
|
kpeter@814
|
361 |
|
kpeter@814
|
362 |
/// \brief Return the flow on the given arc.
|
kpeter@814
|
363 |
///
|
kpeter@814
|
364 |
/// Return the flow on the given arc.
|
kpeter@814
|
365 |
///
|
kpeter@814
|
366 |
/// \pre \ref run() must be called before using this function.
|
kpeter@814
|
367 |
Capacity flow(const Arc& arc) const {
|
kpeter@814
|
368 |
return (*_flow)[arc];
|
kpeter@814
|
369 |
}
|
kpeter@814
|
370 |
|
kpeter@814
|
371 |
/// \brief Return the potential of the given node.
|
kpeter@814
|
372 |
///
|
kpeter@814
|
373 |
/// Return the potential of the given node.
|
kpeter@814
|
374 |
///
|
kpeter@814
|
375 |
/// \pre \ref run() must be called before using this function.
|
kpeter@814
|
376 |
Cost potential(const Node& node) const {
|
kpeter@814
|
377 |
return (*_potential)[node];
|
kpeter@814
|
378 |
}
|
kpeter@814
|
379 |
|
kpeter@814
|
380 |
/// \brief Return the total cost of the found flow.
|
kpeter@814
|
381 |
///
|
kpeter@814
|
382 |
/// Return the total cost of the found flow. The complexity of the
|
kpeter@814
|
383 |
/// function is \f$ O(e) \f$.
|
kpeter@814
|
384 |
///
|
kpeter@814
|
385 |
/// \pre \ref run() must be called before using this function.
|
kpeter@814
|
386 |
Cost totalCost() const {
|
kpeter@814
|
387 |
Cost c = 0;
|
kpeter@814
|
388 |
for (ArcIt e(_graph); e != INVALID; ++e)
|
kpeter@814
|
389 |
c += (*_flow)[e] * _cost[e];
|
kpeter@814
|
390 |
return c;
|
kpeter@814
|
391 |
}
|
kpeter@814
|
392 |
|
kpeter@814
|
393 |
/// @}
|
kpeter@814
|
394 |
|
kpeter@814
|
395 |
private:
|
kpeter@814
|
396 |
|
kpeter@814
|
397 |
/// Initialize the algorithm.
|
kpeter@814
|
398 |
bool init() {
|
kpeter@814
|
399 |
if (!_valid_supply) return false;
|
kpeter@814
|
400 |
|
kpeter@814
|
401 |
// Initializing flow and potential maps
|
kpeter@814
|
402 |
if (!_flow) {
|
kpeter@814
|
403 |
_flow = new FlowMap(_graph);
|
kpeter@814
|
404 |
_local_flow = true;
|
kpeter@814
|
405 |
}
|
kpeter@814
|
406 |
if (!_potential) {
|
kpeter@814
|
407 |
_potential = new PotentialMap(_graph);
|
kpeter@814
|
408 |
_local_potential = true;
|
kpeter@814
|
409 |
}
|
kpeter@814
|
410 |
|
kpeter@814
|
411 |
_res_graph = new ResDigraph(_graph, _capacity, *_flow);
|
kpeter@814
|
412 |
|
kpeter@814
|
413 |
// Finding a feasible flow using Circulation
|
kpeter@814
|
414 |
Circulation< Digraph, ConstMap<Arc, Capacity>, CapacityArcMap,
|
kpeter@814
|
415 |
SupplyMap >
|
kpeter@814
|
416 |
circulation( _graph, constMap<Arc>(Capacity(0)), _capacity,
|
kpeter@814
|
417 |
_supply );
|
kpeter@814
|
418 |
return circulation.flowMap(*_flow).run();
|
kpeter@814
|
419 |
}
|
kpeter@814
|
420 |
|
kpeter@814
|
421 |
bool start(bool min_mean_cc) {
|
kpeter@814
|
422 |
if (min_mean_cc)
|
kpeter@814
|
423 |
startMinMean();
|
kpeter@814
|
424 |
else
|
kpeter@814
|
425 |
start();
|
kpeter@814
|
426 |
|
kpeter@814
|
427 |
// Handling non-zero lower bounds
|
kpeter@814
|
428 |
if (_lower) {
|
kpeter@814
|
429 |
for (ArcIt e(_graph); e != INVALID; ++e)
|
kpeter@814
|
430 |
(*_flow)[e] += (*_lower)[e];
|
kpeter@814
|
431 |
}
|
kpeter@814
|
432 |
return true;
|
kpeter@814
|
433 |
}
|
kpeter@814
|
434 |
|
kpeter@814
|
435 |
/// \brief Execute the algorithm using \ref BellmanFord.
|
kpeter@814
|
436 |
///
|
kpeter@814
|
437 |
/// Execute the algorithm using the \ref BellmanFord
|
kpeter@814
|
438 |
/// "Bellman-Ford" algorithm for negative cycle detection with
|
kpeter@814
|
439 |
/// successively larger limit for the number of iterations.
|
kpeter@814
|
440 |
void start() {
|
kpeter@814
|
441 |
typename BellmanFord<ResDigraph, ResidualCostMap>::PredMap pred(*_res_graph);
|
kpeter@814
|
442 |
typename ResDigraph::template NodeMap<int> visited(*_res_graph);
|
kpeter@814
|
443 |
std::vector<ResArc> cycle;
|
kpeter@814
|
444 |
int node_num = countNodes(_graph);
|
kpeter@814
|
445 |
|
kpeter@814
|
446 |
int length_bound = BF_FIRST_LIMIT;
|
kpeter@814
|
447 |
bool optimal = false;
|
kpeter@814
|
448 |
while (!optimal) {
|
kpeter@814
|
449 |
BellmanFord<ResDigraph, ResidualCostMap> bf(*_res_graph, _res_cost);
|
kpeter@814
|
450 |
bf.predMap(pred);
|
kpeter@814
|
451 |
bf.init(0);
|
kpeter@814
|
452 |
int iter_num = 0;
|
kpeter@814
|
453 |
bool cycle_found = false;
|
kpeter@814
|
454 |
while (!cycle_found) {
|
kpeter@814
|
455 |
int curr_iter_num = iter_num + length_bound <= node_num ?
|
kpeter@814
|
456 |
length_bound : node_num - iter_num;
|
kpeter@814
|
457 |
iter_num += curr_iter_num;
|
kpeter@814
|
458 |
int real_iter_num = curr_iter_num;
|
kpeter@814
|
459 |
for (int i = 0; i < curr_iter_num; ++i) {
|
kpeter@814
|
460 |
if (bf.processNextWeakRound()) {
|
kpeter@814
|
461 |
real_iter_num = i;
|
kpeter@814
|
462 |
break;
|
kpeter@814
|
463 |
}
|
kpeter@814
|
464 |
}
|
kpeter@814
|
465 |
if (real_iter_num < curr_iter_num) {
|
kpeter@814
|
466 |
// Optimal flow is found
|
kpeter@814
|
467 |
optimal = true;
|
kpeter@814
|
468 |
// Setting node potentials
|
kpeter@814
|
469 |
for (NodeIt n(_graph); n != INVALID; ++n)
|
kpeter@814
|
470 |
(*_potential)[n] = bf.dist(n);
|
kpeter@814
|
471 |
break;
|
kpeter@814
|
472 |
} else {
|
kpeter@814
|
473 |
// Searching for node disjoint negative cycles
|
kpeter@814
|
474 |
for (ResNodeIt n(*_res_graph); n != INVALID; ++n)
|
kpeter@814
|
475 |
visited[n] = 0;
|
kpeter@814
|
476 |
int id = 0;
|
kpeter@814
|
477 |
for (ResNodeIt n(*_res_graph); n != INVALID; ++n) {
|
kpeter@814
|
478 |
if (visited[n] > 0) continue;
|
kpeter@814
|
479 |
visited[n] = ++id;
|
kpeter@814
|
480 |
ResNode u = pred[n] == INVALID ?
|
kpeter@814
|
481 |
INVALID : _res_graph->source(pred[n]);
|
kpeter@814
|
482 |
while (u != INVALID && visited[u] == 0) {
|
kpeter@814
|
483 |
visited[u] = id;
|
kpeter@814
|
484 |
u = pred[u] == INVALID ?
|
kpeter@814
|
485 |
INVALID : _res_graph->source(pred[u]);
|
kpeter@814
|
486 |
}
|
kpeter@814
|
487 |
if (u != INVALID && visited[u] == id) {
|
kpeter@814
|
488 |
// Finding the negative cycle
|
kpeter@814
|
489 |
cycle_found = true;
|
kpeter@814
|
490 |
cycle.clear();
|
kpeter@814
|
491 |
ResArc e = pred[u];
|
kpeter@814
|
492 |
cycle.push_back(e);
|
kpeter@814
|
493 |
Capacity d = _res_graph->residualCapacity(e);
|
kpeter@814
|
494 |
while (_res_graph->source(e) != u) {
|
kpeter@814
|
495 |
cycle.push_back(e = pred[_res_graph->source(e)]);
|
kpeter@814
|
496 |
if (_res_graph->residualCapacity(e) < d)
|
kpeter@814
|
497 |
d = _res_graph->residualCapacity(e);
|
kpeter@814
|
498 |
}
|
kpeter@814
|
499 |
|
kpeter@814
|
500 |
// Augmenting along the cycle
|
kpeter@814
|
501 |
for (int i = 0; i < int(cycle.size()); ++i)
|
kpeter@814
|
502 |
_res_graph->augment(cycle[i], d);
|
kpeter@814
|
503 |
}
|
kpeter@814
|
504 |
}
|
kpeter@814
|
505 |
}
|
kpeter@814
|
506 |
|
kpeter@814
|
507 |
if (!cycle_found)
|
kpeter@814
|
508 |
length_bound = length_bound * BF_LIMIT_FACTOR / 100;
|
kpeter@814
|
509 |
}
|
kpeter@814
|
510 |
}
|
kpeter@814
|
511 |
}
|
kpeter@814
|
512 |
|
kpeter@814
|
513 |
/// \brief Execute the algorithm using \ref Howard.
|
kpeter@814
|
514 |
///
|
kpeter@814
|
515 |
/// Execute the algorithm using \ref Howard for negative
|
kpeter@814
|
516 |
/// cycle detection.
|
kpeter@814
|
517 |
void startMinMean() {
|
kpeter@814
|
518 |
typedef Path<ResDigraph> ResPath;
|
kpeter@814
|
519 |
Howard<ResDigraph, ResidualCostMap> mmc(*_res_graph, _res_cost);
|
kpeter@814
|
520 |
ResPath cycle;
|
kpeter@814
|
521 |
|
kpeter@814
|
522 |
mmc.cycle(cycle);
|
kpeter@814
|
523 |
if (mmc.findMinMean()) {
|
kpeter@814
|
524 |
while (mmc.cycleLength() < 0) {
|
kpeter@814
|
525 |
// Finding the cycle
|
kpeter@814
|
526 |
mmc.findCycle();
|
kpeter@814
|
527 |
|
kpeter@814
|
528 |
// Finding the largest flow amount that can be augmented
|
kpeter@814
|
529 |
// along the cycle
|
kpeter@814
|
530 |
Capacity delta = 0;
|
kpeter@814
|
531 |
for (typename ResPath::ArcIt e(cycle); e != INVALID; ++e) {
|
kpeter@814
|
532 |
if (delta == 0 || _res_graph->residualCapacity(e) < delta)
|
kpeter@814
|
533 |
delta = _res_graph->residualCapacity(e);
|
kpeter@814
|
534 |
}
|
kpeter@814
|
535 |
|
kpeter@814
|
536 |
// Augmenting along the cycle
|
kpeter@814
|
537 |
for (typename ResPath::ArcIt e(cycle); e != INVALID; ++e)
|
kpeter@814
|
538 |
_res_graph->augment(e, delta);
|
kpeter@814
|
539 |
|
kpeter@814
|
540 |
// Finding the minimum cycle mean for the modified residual
|
kpeter@814
|
541 |
// digraph
|
kpeter@814
|
542 |
if (!mmc.findMinMean()) break;
|
kpeter@814
|
543 |
}
|
kpeter@814
|
544 |
}
|
kpeter@814
|
545 |
|
kpeter@814
|
546 |
// Computing node potentials
|
kpeter@814
|
547 |
BellmanFord<ResDigraph, ResidualCostMap> bf(*_res_graph, _res_cost);
|
kpeter@814
|
548 |
bf.init(0); bf.start();
|
kpeter@814
|
549 |
for (NodeIt n(_graph); n != INVALID; ++n)
|
kpeter@814
|
550 |
(*_potential)[n] = bf.dist(n);
|
kpeter@814
|
551 |
}
|
kpeter@814
|
552 |
|
kpeter@814
|
553 |
}; //class CycleCanceling
|
kpeter@814
|
554 |
|
kpeter@814
|
555 |
///@}
|
kpeter@814
|
556 |
|
kpeter@814
|
557 |
} //namespace lemon
|
kpeter@814
|
558 |
|
kpeter@814
|
559 |
#endif //LEMON_CYCLE_CANCELING_H
|