0
3
2
1 |
/* -*- mode: C++; indent-tabs-mode: nil; -*- |
|
2 |
* |
|
3 |
* This file is a part of LEMON, a generic C++ optimization library. |
|
4 |
* |
|
5 |
* Copyright (C) 2003-2008 |
|
6 |
* Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport |
|
7 |
* (Egervary Research Group on Combinatorial Optimization, EGRES). |
|
8 |
* |
|
9 |
* Permission to use, modify and distribute this software is granted |
|
10 |
* provided that this copyright notice appears in all copies. For |
|
11 |
* precise terms see the accompanying LICENSE file. |
|
12 |
* |
|
13 |
* This software is provided "AS IS" with no warranty of any kind, |
|
14 |
* express or implied, and with no claim as to its suitability for any |
|
15 |
* purpose. |
|
16 |
* |
|
17 |
*/ |
|
18 |
|
|
19 |
#ifndef LEMON_CIRCULATION_H |
|
20 |
#define LEMON_CIRCULATION_H |
|
21 |
|
|
22 |
#include <lemon/tolerance.h> |
|
23 |
#include <lemon/elevator.h> |
|
24 |
|
|
25 |
///\ingroup max_flow |
|
26 |
///\file |
|
27 |
///\brief Push-relabel algorithm for finding a feasible circulation. |
|
28 |
/// |
|
29 |
namespace lemon { |
|
30 |
|
|
31 |
/// \brief Default traits class of Circulation class. |
|
32 |
/// |
|
33 |
/// Default traits class of Circulation class. |
|
34 |
/// \tparam _Diraph Digraph type. |
|
35 |
/// \tparam _LCapMap Lower bound capacity map type. |
|
36 |
/// \tparam _UCapMap Upper bound capacity map type. |
|
37 |
/// \tparam _DeltaMap Delta map type. |
|
38 |
template <typename _Diraph, typename _LCapMap, |
|
39 |
typename _UCapMap, typename _DeltaMap> |
|
40 |
struct CirculationDefaultTraits { |
|
41 |
|
|
42 |
/// \brief The type of the digraph the algorithm runs on. |
|
43 |
typedef _Diraph Digraph; |
|
44 |
|
|
45 |
/// \brief The type of the map that stores the circulation lower |
|
46 |
/// bound. |
|
47 |
/// |
|
48 |
/// The type of the map that stores the circulation lower bound. |
|
49 |
/// It must meet the \ref concepts::ReadMap "ReadMap" concept. |
|
50 |
typedef _LCapMap LCapMap; |
|
51 |
|
|
52 |
/// \brief The type of the map that stores the circulation upper |
|
53 |
/// bound. |
|
54 |
/// |
|
55 |
/// The type of the map that stores the circulation upper bound. |
|
56 |
/// It must meet the \ref concepts::ReadMap "ReadMap" concept. |
|
57 |
typedef _UCapMap UCapMap; |
|
58 |
|
|
59 |
/// \brief The type of the map that stores the lower bound for |
|
60 |
/// the supply of the nodes. |
|
61 |
/// |
|
62 |
/// The type of the map that stores the lower bound for the supply |
|
63 |
/// of the nodes. It must meet the \ref concepts::ReadMap "ReadMap" |
|
64 |
/// concept. |
|
65 |
typedef _DeltaMap DeltaMap; |
|
66 |
|
|
67 |
/// \brief The type of the flow values. |
|
68 |
typedef typename DeltaMap::Value Value; |
|
69 |
|
|
70 |
/// \brief The type of the map that stores the flow values. |
|
71 |
/// |
|
72 |
/// The type of the map that stores the flow values. |
|
73 |
/// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept. |
|
74 |
typedef typename Digraph::template ArcMap<Value> FlowMap; |
|
75 |
|
|
76 |
/// \brief Instantiates a FlowMap. |
|
77 |
/// |
|
78 |
/// This function instantiates a \ref FlowMap. |
|
79 |
/// \param digraph The digraph, to which we would like to define |
|
80 |
/// the flow map. |
|
81 |
static FlowMap* createFlowMap(const Digraph& digraph) { |
|
82 |
return new FlowMap(digraph); |
|
83 |
} |
|
84 |
|
|
85 |
/// \brief The elevator type used by the algorithm. |
|
86 |
/// |
|
87 |
/// The elevator type used by the algorithm. |
|
88 |
/// |
|
89 |
/// \sa Elevator |
|
90 |
/// \sa LinkedElevator |
|
91 |
typedef lemon::Elevator<Digraph, typename Digraph::Node> Elevator; |
|
92 |
|
|
93 |
/// \brief Instantiates an Elevator. |
|
94 |
/// |
|
95 |
/// This function instantiates an \ref Elevator. |
|
96 |
/// \param digraph The digraph, to which we would like to define |
|
97 |
/// the elevator. |
|
98 |
/// \param max_level The maximum level of the elevator. |
|
99 |
static Elevator* createElevator(const Digraph& digraph, int max_level) { |
|
100 |
return new Elevator(digraph, max_level); |
|
101 |
} |
|
102 |
|
|
103 |
/// \brief The tolerance used by the algorithm |
|
104 |
/// |
|
105 |
/// The tolerance used by the algorithm to handle inexact computation. |
|
106 |
typedef lemon::Tolerance<Value> Tolerance; |
|
107 |
|
|
108 |
}; |
|
109 |
|
|
110 |
/** |
|
111 |
\brief Push-relabel algorithm for the network circulation problem. |
|
112 |
|
|
113 |
\ingroup max_flow |
|
114 |
This class implements a push-relabel algorithm for the network |
|
115 |
circulation problem. |
|
116 |
It is to find a feasible circulation when lower and upper bounds |
|
117 |
are given for the flow values on the arcs and lower bounds |
|
118 |
are given for the supply values of the nodes. |
|
119 |
|
|
120 |
The exact formulation of this problem is the following. |
|
121 |
Let \f$G=(V,A)\f$ be a digraph, |
|
122 |
\f$lower, upper: A\rightarrow\mathbf{R}^+_0\f$, |
|
123 |
\f$delta: V\rightarrow\mathbf{R}\f$. Find a feasible circulation |
|
124 |
\f$f: A\rightarrow\mathbf{R}^+_0\f$ so that |
|
125 |
\f[ \sum_{a\in\delta_{out}(v)} f(a) - \sum_{a\in\delta_{in}(v)} f(a) |
|
126 |
\geq delta(v) \quad \forall v\in V, \f] |
|
127 |
\f[ lower(a)\leq f(a) \leq upper(a) \quad \forall a\in A. \f] |
|
128 |
\note \f$delta(v)\f$ specifies a lower bound for the supply of node |
|
129 |
\f$v\f$. It can be either positive or negative, however note that |
|
130 |
\f$\sum_{v\in V}delta(v)\f$ should be zero or negative in order to |
|
131 |
have a feasible solution. |
|
132 |
|
|
133 |
\note A special case of this problem is when |
|
134 |
\f$\sum_{v\in V}delta(v) = 0\f$. Then the supply of each node \f$v\f$ |
|
135 |
will be \e equal \e to \f$delta(v)\f$, if a circulation can be found. |
|
136 |
Thus a feasible solution for the |
|
137 |
\ref min_cost_flow "minimum cost flow" problem can be calculated |
|
138 |
in this way. |
|
139 |
|
|
140 |
\tparam _Digraph The type of the digraph the algorithm runs on. |
|
141 |
\tparam _LCapMap The type of the lower bound capacity map. The default |
|
142 |
map type is \ref concepts::Digraph::ArcMap "_Digraph::ArcMap<int>". |
|
143 |
\tparam _UCapMap The type of the upper bound capacity map. The default |
|
144 |
map type is \c _LCapMap. |
|
145 |
\tparam _DeltaMap The type of the map that stores the lower bound |
|
146 |
for the supply of the nodes. The default map type is |
|
147 |
\c _Digraph::ArcMap<_UCapMap::Value>. |
|
148 |
*/ |
|
149 |
#ifdef DOXYGEN |
|
150 |
template< typename _Digraph, |
|
151 |
typename _LCapMap, |
|
152 |
typename _UCapMap, |
|
153 |
typename _DeltaMap, |
|
154 |
typename _Traits > |
|
155 |
#else |
|
156 |
template< typename _Digraph, |
|
157 |
typename _LCapMap = typename _Digraph::template ArcMap<int>, |
|
158 |
typename _UCapMap = _LCapMap, |
|
159 |
typename _DeltaMap = typename _Digraph:: |
|
160 |
template NodeMap<typename _UCapMap::Value>, |
|
161 |
typename _Traits=CirculationDefaultTraits<_Digraph, _LCapMap, |
|
162 |
_UCapMap, _DeltaMap> > |
|
163 |
#endif |
|
164 |
class Circulation { |
|
165 |
public: |
|
166 |
|
|
167 |
///The \ref CirculationDefaultTraits "traits class" of the algorithm. |
|
168 |
typedef _Traits Traits; |
|
169 |
///The type of the digraph the algorithm runs on. |
|
170 |
typedef typename Traits::Digraph Digraph; |
|
171 |
///The type of the flow values. |
|
172 |
typedef typename Traits::Value Value; |
|
173 |
|
|
174 |
/// The type of the lower bound capacity map. |
|
175 |
typedef typename Traits::LCapMap LCapMap; |
|
176 |
/// The type of the upper bound capacity map. |
|
177 |
typedef typename Traits::UCapMap UCapMap; |
|
178 |
/// \brief The type of the map that stores the lower bound for |
|
179 |
/// the supply of the nodes. |
|
180 |
typedef typename Traits::DeltaMap DeltaMap; |
|
181 |
///The type of the flow map. |
|
182 |
typedef typename Traits::FlowMap FlowMap; |
|
183 |
|
|
184 |
///The type of the elevator. |
|
185 |
typedef typename Traits::Elevator Elevator; |
|
186 |
///The type of the tolerance. |
|
187 |
typedef typename Traits::Tolerance Tolerance; |
|
188 |
|
|
189 |
private: |
|
190 |
|
|
191 |
TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); |
|
192 |
|
|
193 |
const Digraph &_g; |
|
194 |
int _node_num; |
|
195 |
|
|
196 |
const LCapMap *_lo; |
|
197 |
const UCapMap *_up; |
|
198 |
const DeltaMap *_delta; |
|
199 |
|
|
200 |
FlowMap *_flow; |
|
201 |
bool _local_flow; |
|
202 |
|
|
203 |
Elevator* _level; |
|
204 |
bool _local_level; |
|
205 |
|
|
206 |
typedef typename Digraph::template NodeMap<Value> ExcessMap; |
|
207 |
ExcessMap* _excess; |
|
208 |
|
|
209 |
Tolerance _tol; |
|
210 |
int _el; |
|
211 |
|
|
212 |
public: |
|
213 |
|
|
214 |
typedef Circulation Create; |
|
215 |
|
|
216 |
///\name Named Template Parameters |
|
217 |
|
|
218 |
///@{ |
|
219 |
|
|
220 |
template <typename _FlowMap> |
|
221 |
struct SetFlowMapTraits : public Traits { |
|
222 |
typedef _FlowMap FlowMap; |
|
223 |
static FlowMap *createFlowMap(const Digraph&) { |
|
224 |
LEMON_ASSERT(false, "FlowMap is not initialized"); |
|
225 |
return 0; // ignore warnings |
|
226 |
} |
|
227 |
}; |
|
228 |
|
|
229 |
/// \brief \ref named-templ-param "Named parameter" for setting |
|
230 |
/// FlowMap type |
|
231 |
/// |
|
232 |
/// \ref named-templ-param "Named parameter" for setting FlowMap |
|
233 |
/// type. |
|
234 |
template <typename _FlowMap> |
|
235 |
struct SetFlowMap |
|
236 |
: public Circulation<Digraph, LCapMap, UCapMap, DeltaMap, |
|
237 |
SetFlowMapTraits<_FlowMap> > { |
|
238 |
typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap, |
|
239 |
SetFlowMapTraits<_FlowMap> > Create; |
|
240 |
}; |
|
241 |
|
|
242 |
template <typename _Elevator> |
|
243 |
struct SetElevatorTraits : public Traits { |
|
244 |
typedef _Elevator Elevator; |
|
245 |
static Elevator *createElevator(const Digraph&, int) { |
|
246 |
LEMON_ASSERT(false, "Elevator is not initialized"); |
|
247 |
return 0; // ignore warnings |
|
248 |
} |
|
249 |
}; |
|
250 |
|
|
251 |
/// \brief \ref named-templ-param "Named parameter" for setting |
|
252 |
/// Elevator type |
|
253 |
/// |
|
254 |
/// \ref named-templ-param "Named parameter" for setting Elevator |
|
255 |
/// type. If this named parameter is used, then an external |
|
256 |
/// elevator object must be passed to the algorithm using the |
|
257 |
/// \ref elevator(Elevator&) "elevator()" function before calling |
|
258 |
/// \ref run() or \ref init(). |
|
259 |
/// \sa SetStandardElevator |
|
260 |
template <typename _Elevator> |
|
261 |
struct SetElevator |
|
262 |
: public Circulation<Digraph, LCapMap, UCapMap, DeltaMap, |
|
263 |
SetElevatorTraits<_Elevator> > { |
|
264 |
typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap, |
|
265 |
SetElevatorTraits<_Elevator> > Create; |
|
266 |
}; |
|
267 |
|
|
268 |
template <typename _Elevator> |
|
269 |
struct SetStandardElevatorTraits : public Traits { |
|
270 |
typedef _Elevator Elevator; |
|
271 |
static Elevator *createElevator(const Digraph& digraph, int max_level) { |
|
272 |
return new Elevator(digraph, max_level); |
|
273 |
} |
|
274 |
}; |
|
275 |
|
|
276 |
/// \brief \ref named-templ-param "Named parameter" for setting |
|
277 |
/// Elevator type with automatic allocation |
|
278 |
/// |
|
279 |
/// \ref named-templ-param "Named parameter" for setting Elevator |
|
280 |
/// type with automatic allocation. |
|
281 |
/// The Elevator should have standard constructor interface to be |
|
282 |
/// able to automatically created by the algorithm (i.e. the |
|
283 |
/// digraph and the maximum level should be passed to it). |
|
284 |
/// However an external elevator object could also be passed to the |
|
285 |
/// algorithm with the \ref elevator(Elevator&) "elevator()" function |
|
286 |
/// before calling \ref run() or \ref init(). |
|
287 |
/// \sa SetElevator |
|
288 |
template <typename _Elevator> |
|
289 |
struct SetStandardElevator |
|
290 |
: public Circulation<Digraph, LCapMap, UCapMap, DeltaMap, |
|
291 |
SetStandardElevatorTraits<_Elevator> > { |
|
292 |
typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap, |
|
293 |
SetStandardElevatorTraits<_Elevator> > Create; |
|
294 |
}; |
|
295 |
|
|
296 |
/// @} |
|
297 |
|
|
298 |
protected: |
|
299 |
|
|
300 |
Circulation() {} |
|
301 |
|
|
302 |
public: |
|
303 |
|
|
304 |
/// The constructor of the class. |
|
305 |
|
|
306 |
/// The constructor of the class. |
|
307 |
/// \param g The digraph the algorithm runs on. |
|
308 |
/// \param lo The lower bound capacity of the arcs. |
|
309 |
/// \param up The upper bound capacity of the arcs. |
|
310 |
/// \param delta The lower bound for the supply of the nodes. |
|
311 |
Circulation(const Digraph &g,const LCapMap &lo, |
|
312 |
const UCapMap &up,const DeltaMap &delta) |
|
313 |
: _g(g), _node_num(), |
|
314 |
_lo(&lo),_up(&up),_delta(&delta),_flow(0),_local_flow(false), |
|
315 |
_level(0), _local_level(false), _excess(0), _el() {} |
|
316 |
|
|
317 |
/// Destructor. |
|
318 |
~Circulation() { |
|
319 |
destroyStructures(); |
|
320 |
} |
|
321 |
|
|
322 |
|
|
323 |
private: |
|
324 |
|
|
325 |
void createStructures() { |
|
326 |
_node_num = _el = countNodes(_g); |
|
327 |
|
|
328 |
if (!_flow) { |
|
329 |
_flow = Traits::createFlowMap(_g); |
|
330 |
_local_flow = true; |
|
331 |
} |
|
332 |
if (!_level) { |
|
333 |
_level = Traits::createElevator(_g, _node_num); |
|
334 |
_local_level = true; |
|
335 |
} |
|
336 |
if (!_excess) { |
|
337 |
_excess = new ExcessMap(_g); |
|
338 |
} |
|
339 |
} |
|
340 |
|
|
341 |
void destroyStructures() { |
|
342 |
if (_local_flow) { |
|
343 |
delete _flow; |
|
344 |
} |
|
345 |
if (_local_level) { |
|
346 |
delete _level; |
|
347 |
} |
|
348 |
if (_excess) { |
|
349 |
delete _excess; |
|
350 |
} |
|
351 |
} |
|
352 |
|
|
353 |
public: |
|
354 |
|
|
355 |
/// Sets the lower bound capacity map. |
|
356 |
|
|
357 |
/// Sets the lower bound capacity map. |
|
358 |
/// \return <tt>(*this)</tt> |
|
359 |
Circulation& lowerCapMap(const LCapMap& map) { |
|
360 |
_lo = ↦ |
|
361 |
return *this; |
|
362 |
} |
|
363 |
|
|
364 |
/// Sets the upper bound capacity map. |
|
365 |
|
|
366 |
/// Sets the upper bound capacity map. |
|
367 |
/// \return <tt>(*this)</tt> |
|
368 |
Circulation& upperCapMap(const LCapMap& map) { |
|
369 |
_up = ↦ |
|
370 |
return *this; |
|
371 |
} |
|
372 |
|
|
373 |
/// Sets the lower bound map for the supply of the nodes. |
|
374 |
|
|
375 |
/// Sets the lower bound map for the supply of the nodes. |
|
376 |
/// \return <tt>(*this)</tt> |
|
377 |
Circulation& deltaMap(const DeltaMap& map) { |
|
378 |
_delta = ↦ |
|
379 |
return *this; |
|
380 |
} |
|
381 |
|
|
382 |
/// \brief Sets the flow map. |
|
383 |
/// |
|
384 |
/// Sets the flow map. |
|
385 |
/// If you don't use this function before calling \ref run() or |
|
386 |
/// \ref init(), an instance will be allocated automatically. |
|
387 |
/// The destructor deallocates this automatically allocated map, |
|
388 |
/// of course. |
|
389 |
/// \return <tt>(*this)</tt> |
|
390 |
Circulation& flowMap(FlowMap& map) { |
|
391 |
if (_local_flow) { |
|
392 |
delete _flow; |
|
393 |
_local_flow = false; |
|
394 |
} |
|
395 |
_flow = ↦ |
|
396 |
return *this; |
|
397 |
} |
|
398 |
|
|
399 |
/// \brief Sets the elevator used by algorithm. |
|
400 |
/// |
|
401 |
/// Sets the elevator used by algorithm. |
|
402 |
/// If you don't use this function before calling \ref run() or |
|
403 |
/// \ref init(), an instance will be allocated automatically. |
|
404 |
/// The destructor deallocates this automatically allocated elevator, |
|
405 |
/// of course. |
|
406 |
/// \return <tt>(*this)</tt> |
|
407 |
Circulation& elevator(Elevator& elevator) { |
|
408 |
if (_local_level) { |
|
409 |
delete _level; |
|
410 |
_local_level = false; |
|
411 |
} |
|
412 |
_level = &elevator; |
|
413 |
return *this; |
|
414 |
} |
|
415 |
|
|
416 |
/// \brief Returns a const reference to the elevator. |
|
417 |
/// |
|
418 |
/// Returns a const reference to the elevator. |
|
419 |
/// |
|
420 |
/// \pre Either \ref run() or \ref init() must be called before |
|
421 |
/// using this function. |
|
422 |
const Elevator& elevator() { |
|
423 |
return *_level; |
|
424 |
} |
|
425 |
|
|
426 |
/// \brief Sets the tolerance used by algorithm. |
|
427 |
/// |
|
428 |
/// Sets the tolerance used by algorithm. |
|
429 |
Circulation& tolerance(const Tolerance& tolerance) const { |
|
430 |
_tol = tolerance; |
|
431 |
return *this; |
|
432 |
} |
|
433 |
|
|
434 |
/// \brief Returns a const reference to the tolerance. |
|
435 |
/// |
|
436 |
/// Returns a const reference to the tolerance. |
|
437 |
const Tolerance& tolerance() const { |
|
438 |
return tolerance; |
|
439 |
} |
|
440 |
|
|
441 |
/// \name Execution Control |
|
442 |
/// The simplest way to execute the algorithm is to call \ref run().\n |
|
443 |
/// If you need more control on the initial solution or the execution, |
|
444 |
/// first you have to call one of the \ref init() functions, then |
|
445 |
/// the \ref start() function. |
|
446 |
|
|
447 |
///@{ |
|
448 |
|
|
449 |
/// Initializes the internal data structures. |
|
450 |
|
|
451 |
/// Initializes the internal data structures and sets all flow values |
|
452 |
/// to the lower bound. |
|
453 |
void init() |
|
454 |
{ |
|
455 |
createStructures(); |
|
456 |
|
|
457 |
for(NodeIt n(_g);n!=INVALID;++n) { |
|
458 |
_excess->set(n, (*_delta)[n]); |
|
459 |
} |
|
460 |
|
|
461 |
for (ArcIt e(_g);e!=INVALID;++e) { |
|
462 |
_flow->set(e, (*_lo)[e]); |
|
463 |
_excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_flow)[e]); |
|
464 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_flow)[e]); |
|
465 |
} |
|
466 |
|
|
467 |
// global relabeling tested, but in general case it provides |
|
468 |
// worse performance for random digraphs |
|
469 |
_level->initStart(); |
|
470 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
471 |
_level->initAddItem(n); |
|
472 |
_level->initFinish(); |
|
473 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
474 |
if(_tol.positive((*_excess)[n])) |
|
475 |
_level->activate(n); |
|
476 |
} |
|
477 |
|
|
478 |
/// Initializes the internal data structures using a greedy approach. |
|
479 |
|
|
480 |
/// Initializes the internal data structures using a greedy approach |
|
481 |
/// to construct the initial solution. |
|
482 |
void greedyInit() |
|
483 |
{ |
|
484 |
createStructures(); |
|
485 |
|
|
486 |
for(NodeIt n(_g);n!=INVALID;++n) { |
|
487 |
_excess->set(n, (*_delta)[n]); |
|
488 |
} |
|
489 |
|
|
490 |
for (ArcIt e(_g);e!=INVALID;++e) { |
|
491 |
if (!_tol.positive((*_excess)[_g.target(e)] + (*_up)[e])) { |
|
492 |
_flow->set(e, (*_up)[e]); |
|
493 |
_excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_up)[e]); |
|
494 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_up)[e]); |
|
495 |
} else if (_tol.positive((*_excess)[_g.target(e)] + (*_lo)[e])) { |
|
496 |
_flow->set(e, (*_lo)[e]); |
|
497 |
_excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_lo)[e]); |
|
498 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_lo)[e]); |
|
499 |
} else { |
|
500 |
Value fc = -(*_excess)[_g.target(e)]; |
|
501 |
_flow->set(e, fc); |
|
502 |
_excess->set(_g.target(e), 0); |
|
503 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - fc); |
|
504 |
} |
|
505 |
} |
|
506 |
|
|
507 |
_level->initStart(); |
|
508 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
509 |
_level->initAddItem(n); |
|
510 |
_level->initFinish(); |
|
511 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
512 |
if(_tol.positive((*_excess)[n])) |
|
513 |
_level->activate(n); |
|
514 |
} |
|
515 |
|
|
516 |
///Executes the algorithm |
|
517 |
|
|
518 |
///This function executes the algorithm. |
|
519 |
/// |
|
520 |
///\return \c true if a feasible circulation is found. |
|
521 |
/// |
|
522 |
///\sa barrier() |
|
523 |
///\sa barrierMap() |
|
524 |
bool start() |
|
525 |
{ |
|
526 |
|
|
527 |
Node act; |
|
528 |
Node bact=INVALID; |
|
529 |
Node last_activated=INVALID; |
|
530 |
while((act=_level->highestActive())!=INVALID) { |
|
531 |
int actlevel=(*_level)[act]; |
|
532 |
int mlevel=_node_num; |
|
533 |
Value exc=(*_excess)[act]; |
|
534 |
|
|
535 |
for(OutArcIt e(_g,act);e!=INVALID; ++e) { |
|
536 |
Node v = _g.target(e); |
|
537 |
Value fc=(*_up)[e]-(*_flow)[e]; |
|
538 |
if(!_tol.positive(fc)) continue; |
|
539 |
if((*_level)[v]<actlevel) { |
|
540 |
if(!_tol.less(fc, exc)) { |
|
541 |
_flow->set(e, (*_flow)[e] + exc); |
|
542 |
_excess->set(v, (*_excess)[v] + exc); |
|
543 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
|
544 |
_level->activate(v); |
|
545 |
_excess->set(act,0); |
|
546 |
_level->deactivate(act); |
|
547 |
goto next_l; |
|
548 |
} |
|
549 |
else { |
|
550 |
_flow->set(e, (*_up)[e]); |
|
551 |
_excess->set(v, (*_excess)[v] + fc); |
|
552 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
|
553 |
_level->activate(v); |
|
554 |
exc-=fc; |
|
555 |
} |
|
556 |
} |
|
557 |
else if((*_level)[v]<mlevel) mlevel=(*_level)[v]; |
|
558 |
} |
|
559 |
for(InArcIt e(_g,act);e!=INVALID; ++e) { |
|
560 |
Node v = _g.source(e); |
|
561 |
Value fc=(*_flow)[e]-(*_lo)[e]; |
|
562 |
if(!_tol.positive(fc)) continue; |
|
563 |
if((*_level)[v]<actlevel) { |
|
564 |
if(!_tol.less(fc, exc)) { |
|
565 |
_flow->set(e, (*_flow)[e] - exc); |
|
566 |
_excess->set(v, (*_excess)[v] + exc); |
|
567 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
|
568 |
_level->activate(v); |
|
569 |
_excess->set(act,0); |
|
570 |
_level->deactivate(act); |
|
571 |
goto next_l; |
|
572 |
} |
|
573 |
else { |
|
574 |
_flow->set(e, (*_lo)[e]); |
|
575 |
_excess->set(v, (*_excess)[v] + fc); |
|
576 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
|
577 |
_level->activate(v); |
|
578 |
exc-=fc; |
|
579 |
} |
|
580 |
} |
|
581 |
else if((*_level)[v]<mlevel) mlevel=(*_level)[v]; |
|
582 |
} |
|
583 |
|
|
584 |
_excess->set(act, exc); |
|
585 |
if(!_tol.positive(exc)) _level->deactivate(act); |
|
586 |
else if(mlevel==_node_num) { |
|
587 |
_level->liftHighestActiveToTop(); |
|
588 |
_el = _node_num; |
|
589 |
return false; |
|
590 |
} |
|
591 |
else { |
|
592 |
_level->liftHighestActive(mlevel+1); |
|
593 |
if(_level->onLevel(actlevel)==0) { |
|
594 |
_el = actlevel; |
|
595 |
return false; |
|
596 |
} |
|
597 |
} |
|
598 |
next_l: |
|
599 |
; |
|
600 |
} |
|
601 |
return true; |
|
602 |
} |
|
603 |
|
|
604 |
/// Runs the algorithm. |
|
605 |
|
|
606 |
/// This function runs the algorithm. |
|
607 |
/// |
|
608 |
/// \return \c true if a feasible circulation is found. |
|
609 |
/// |
|
610 |
/// \note Apart from the return value, c.run() is just a shortcut of |
|
611 |
/// the following code. |
|
612 |
/// \code |
|
613 |
/// c.greedyInit(); |
|
614 |
/// c.start(); |
|
615 |
/// \endcode |
|
616 |
bool run() { |
|
617 |
greedyInit(); |
|
618 |
return start(); |
|
619 |
} |
|
620 |
|
|
621 |
/// @} |
|
622 |
|
|
623 |
/// \name Query Functions |
|
624 |
/// The results of the circulation algorithm can be obtained using |
|
625 |
/// these functions.\n |
|
626 |
/// Either \ref run() or \ref start() should be called before |
|
627 |
/// using them. |
|
628 |
|
|
629 |
///@{ |
|
630 |
|
|
631 |
/// \brief Returns the flow on the given arc. |
|
632 |
/// |
|
633 |
/// Returns the flow on the given arc. |
|
634 |
/// |
|
635 |
/// \pre Either \ref run() or \ref init() must be called before |
|
636 |
/// using this function. |
|
637 |
Value flow(const Arc& arc) const { |
|
638 |
return (*_flow)[arc]; |
|
639 |
} |
|
640 |
|
|
641 |
/// \brief Returns a const reference to the flow map. |
|
642 |
/// |
|
643 |
/// Returns a const reference to the arc map storing the found flow. |
|
644 |
/// |
|
645 |
/// \pre Either \ref run() or \ref init() must be called before |
|
646 |
/// using this function. |
|
647 |
const FlowMap& flowMap() { |
|
648 |
return *_flow; |
|
649 |
} |
|
650 |
|
|
651 |
/** |
|
652 |
\brief Returns \c true if the given node is in a barrier. |
|
653 |
|
|
654 |
Barrier is a set \e B of nodes for which |
|
655 |
|
|
656 |
\f[ \sum_{a\in\delta_{out}(B)} upper(a) - |
|
657 |
\sum_{a\in\delta_{in}(B)} lower(a) < \sum_{v\in B}delta(v) \f] |
|
658 |
|
|
659 |
holds. The existence of a set with this property prooves that a |
|
660 |
feasible circualtion cannot exist. |
|
661 |
|
|
662 |
This function returns \c true if the given node is in the found |
|
663 |
barrier. If a feasible circulation is found, the function |
|
664 |
gives back \c false for every node. |
|
665 |
|
|
666 |
\pre Either \ref run() or \ref init() must be called before |
|
667 |
using this function. |
|
668 |
|
|
669 |
\sa barrierMap() |
|
670 |
\sa checkBarrier() |
|
671 |
*/ |
|
672 |
bool barrier(const Node& node) |
|
673 |
{ |
|
674 |
return (*_level)[node] >= _el; |
|
675 |
} |
|
676 |
|
|
677 |
/// \brief Gives back a barrier. |
|
678 |
/// |
|
679 |
/// This function sets \c bar to the characteristic vector of the |
|
680 |
/// found barrier. \c bar should be a \ref concepts::WriteMap "writable" |
|
681 |
/// node map with \c bool (or convertible) value type. |
|
682 |
/// |
|
683 |
/// If a feasible circulation is found, the function gives back an |
|
684 |
/// empty set, so \c bar[v] will be \c false for all nodes \c v. |
|
685 |
/// |
|
686 |
/// \note This function calls \ref barrier() for each node, |
|
687 |
/// so it runs in \f$O(n)\f$ time. |
|
688 |
/// |
|
689 |
/// \pre Either \ref run() or \ref init() must be called before |
|
690 |
/// using this function. |
|
691 |
/// |
|
692 |
/// \sa barrier() |
|
693 |
/// \sa checkBarrier() |
|
694 |
template<class BarrierMap> |
|
695 |
void barrierMap(BarrierMap &bar) |
|
696 |
{ |
|
697 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
698 |
bar.set(n, (*_level)[n] >= _el); |
|
699 |
} |
|
700 |
|
|
701 |
/// @} |
|
702 |
|
|
703 |
/// \name Checker Functions |
|
704 |
/// The feasibility of the results can be checked using |
|
705 |
/// these functions.\n |
|
706 |
/// Either \ref run() or \ref start() should be called before |
|
707 |
/// using them. |
|
708 |
|
|
709 |
///@{ |
|
710 |
|
|
711 |
///Check if the found flow is a feasible circulation |
|
712 |
|
|
713 |
///Check if the found flow is a feasible circulation, |
|
714 |
/// |
|
715 |
bool checkFlow() { |
|
716 |
for(ArcIt e(_g);e!=INVALID;++e) |
|
717 |
if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false; |
|
718 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
719 |
{ |
|
720 |
Value dif=-(*_delta)[n]; |
|
721 |
for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e]; |
|
722 |
for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e]; |
|
723 |
if(_tol.negative(dif)) return false; |
|
724 |
} |
|
725 |
return true; |
|
726 |
} |
|
727 |
|
|
728 |
///Check whether or not the last execution provides a barrier |
|
729 |
|
|
730 |
///Check whether or not the last execution provides a barrier. |
|
731 |
///\sa barrier() |
|
732 |
///\sa barrierMap() |
|
733 |
bool checkBarrier() |
|
734 |
{ |
|
735 |
Value delta=0; |
|
736 |
for(NodeIt n(_g);n!=INVALID;++n) |
|
737 |
if(barrier(n)) |
|
738 |
delta-=(*_delta)[n]; |
|
739 |
for(ArcIt e(_g);e!=INVALID;++e) |
|
740 |
{ |
|
741 |
Node s=_g.source(e); |
|
742 |
Node t=_g.target(e); |
|
743 |
if(barrier(s)&&!barrier(t)) delta+=(*_up)[e]; |
|
744 |
else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e]; |
|
745 |
} |
|
746 |
return _tol.negative(delta); |
|
747 |
} |
|
748 |
|
|
749 |
/// @} |
|
750 |
|
|
751 |
}; |
|
752 |
|
|
753 |
} |
|
754 |
|
|
755 |
#endif |
1 |
/* -*- mode: C++; indent-tabs-mode: nil; -*- |
|
2 |
* |
|
3 |
* This file is a part of LEMON, a generic C++ optimization library. |
|
4 |
* |
|
5 |
* Copyright (C) 2003-2008 |
|
6 |
* Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport |
|
7 |
* (Egervary Research Group on Combinatorial Optimization, EGRES). |
|
8 |
* |
|
9 |
* Permission to use, modify and distribute this software is granted |
|
10 |
* provided that this copyright notice appears in all copies. For |
|
11 |
* precise terms see the accompanying LICENSE file. |
|
12 |
* |
|
13 |
* This software is provided "AS IS" with no warranty of any kind, |
|
14 |
* express or implied, and with no claim as to its suitability for any |
|
15 |
* purpose. |
|
16 |
* |
|
17 |
*/ |
|
18 |
|
|
19 |
#include <iostream> |
|
20 |
|
|
21 |
#include "test_tools.h" |
|
22 |
#include <lemon/list_graph.h> |
|
23 |
#include <lemon/circulation.h> |
|
24 |
#include <lemon/lgf_reader.h> |
|
25 |
#include <lemon/concepts/digraph.h> |
|
26 |
#include <lemon/concepts/maps.h> |
|
27 |
|
|
28 |
using namespace lemon; |
|
29 |
|
|
30 |
char test_lgf[] = |
|
31 |
"@nodes\n" |
|
32 |
"label\n" |
|
33 |
"0\n" |
|
34 |
"1\n" |
|
35 |
"2\n" |
|
36 |
"3\n" |
|
37 |
"4\n" |
|
38 |
"5\n" |
|
39 |
"@arcs\n" |
|
40 |
" lcap ucap\n" |
|
41 |
"0 1 2 10\n" |
|
42 |
"0 2 2 6\n" |
|
43 |
"1 3 4 7\n" |
|
44 |
"1 4 0 5\n" |
|
45 |
"2 4 1 3\n" |
|
46 |
"3 5 3 8\n" |
|
47 |
"4 5 3 7\n" |
|
48 |
"@attributes\n" |
|
49 |
"source 0\n" |
|
50 |
"sink 5\n"; |
|
51 |
|
|
52 |
void checkCirculationCompile() |
|
53 |
{ |
|
54 |
typedef int VType; |
|
55 |
typedef concepts::Digraph Digraph; |
|
56 |
|
|
57 |
typedef Digraph::Node Node; |
|
58 |
typedef Digraph::Arc Arc; |
|
59 |
typedef concepts::ReadMap<Arc,VType> CapMap; |
|
60 |
typedef concepts::ReadMap<Node,VType> DeltaMap; |
|
61 |
typedef concepts::ReadWriteMap<Arc,VType> FlowMap; |
|
62 |
typedef concepts::WriteMap<Node,bool> BarrierMap; |
|
63 |
|
|
64 |
typedef Elevator<Digraph, Digraph::Node> Elev; |
|
65 |
typedef LinkedElevator<Digraph, Digraph::Node> LinkedElev; |
|
66 |
|
|
67 |
Digraph g; |
|
68 |
Node n; |
|
69 |
Arc a; |
|
70 |
CapMap lcap, ucap; |
|
71 |
DeltaMap delta; |
|
72 |
FlowMap flow; |
|
73 |
BarrierMap bar; |
|
74 |
|
|
75 |
Circulation<Digraph, CapMap, CapMap, DeltaMap> |
|
76 |
::SetFlowMap<FlowMap> |
|
77 |
::SetElevator<Elev> |
|
78 |
::SetStandardElevator<LinkedElev> |
|
79 |
::Create circ_test(g,lcap,ucap,delta); |
|
80 |
|
|
81 |
circ_test.lowerCapMap(lcap); |
|
82 |
circ_test.upperCapMap(ucap); |
|
83 |
circ_test.deltaMap(delta); |
|
84 |
flow = circ_test.flowMap(); |
|
85 |
circ_test.flowMap(flow); |
|
86 |
|
|
87 |
circ_test.init(); |
|
88 |
circ_test.greedyInit(); |
|
89 |
circ_test.start(); |
|
90 |
circ_test.run(); |
|
91 |
|
|
92 |
circ_test.barrier(n); |
|
93 |
circ_test.barrierMap(bar); |
|
94 |
circ_test.flow(a); |
|
95 |
} |
|
96 |
|
|
97 |
template <class G, class LM, class UM, class DM> |
|
98 |
void checkCirculation(const G& g, const LM& lm, const UM& um, |
|
99 |
const DM& dm, bool find) |
|
100 |
{ |
|
101 |
Circulation<G, LM, UM, DM> circ(g, lm, um, dm); |
|
102 |
bool ret = circ.run(); |
|
103 |
if (find) { |
|
104 |
check(ret, "A feasible solution should have been found."); |
|
105 |
check(circ.checkFlow(), "The found flow is corrupt."); |
|
106 |
check(!circ.checkBarrier(), "A barrier should not have been found."); |
|
107 |
} else { |
|
108 |
check(!ret, "A feasible solution should not have been found."); |
|
109 |
check(circ.checkBarrier(), "The found barrier is corrupt."); |
|
110 |
} |
|
111 |
} |
|
112 |
|
|
113 |
int main (int, char*[]) |
|
114 |
{ |
|
115 |
typedef ListDigraph Digraph; |
|
116 |
DIGRAPH_TYPEDEFS(Digraph); |
|
117 |
|
|
118 |
Digraph g; |
|
119 |
IntArcMap lo(g), up(g); |
|
120 |
IntNodeMap delta(g, 0); |
|
121 |
Node s, t; |
|
122 |
|
|
123 |
std::istringstream input(test_lgf); |
|
124 |
DigraphReader<Digraph>(g,input). |
|
125 |
arcMap("lcap", lo). |
|
126 |
arcMap("ucap", up). |
|
127 |
node("source",s). |
|
128 |
node("sink",t). |
|
129 |
run(); |
|
130 |
|
|
131 |
delta[s] = 7; delta[t] = -7; |
|
132 |
checkCirculation(g, lo, up, delta, true); |
|
133 |
|
|
134 |
delta[s] = 13; delta[t] = -13; |
|
135 |
checkCirculation(g, lo, up, delta, true); |
|
136 |
|
|
137 |
delta[s] = 6; delta[t] = -6; |
|
138 |
checkCirculation(g, lo, up, delta, false); |
|
139 |
|
|
140 |
delta[s] = 14; delta[t] = -14; |
|
141 |
checkCirculation(g, lo, up, delta, false); |
|
142 |
|
|
143 |
delta[s] = 7; delta[t] = -13; |
|
144 |
checkCirculation(g, lo, up, delta, true); |
|
145 |
|
|
146 |
delta[s] = 5; delta[t] = -15; |
|
147 |
checkCirculation(g, lo, up, delta, true); |
|
148 |
|
|
149 |
delta[s] = 10; delta[t] = -11; |
|
150 |
checkCirculation(g, lo, up, delta, true); |
|
151 |
|
|
152 |
delta[s] = 11; delta[t] = -10; |
|
153 |
checkCirculation(g, lo, up, delta, false); |
|
154 |
|
|
155 |
return 0; |
|
156 |
} |
0 comments (0 inline)