0
3
0
1 | 1 |
/* -*- mode: C++; indent-tabs-mode: nil; -*- |
2 | 2 |
* |
3 | 3 |
* This file is a part of LEMON, a generic C++ optimization library. |
4 | 4 |
* |
5 | 5 |
* Copyright (C) 2003-2009 |
6 | 6 |
* Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport |
7 | 7 |
* (Egervary Research Group on Combinatorial Optimization, EGRES). |
8 | 8 |
* |
9 | 9 |
* Permission to use, modify and distribute this software is granted |
10 | 10 |
* provided that this copyright notice appears in all copies. For |
11 | 11 |
* precise terms see the accompanying LICENSE file. |
12 | 12 |
* |
13 | 13 |
* This software is provided "AS IS" with no warranty of any kind, |
14 | 14 |
* express or implied, and with no claim as to its suitability for any |
15 | 15 |
* purpose. |
16 | 16 |
* |
17 | 17 |
*/ |
18 | 18 |
|
19 | 19 |
#ifndef LEMON_CIRCULATION_H |
20 | 20 |
#define LEMON_CIRCULATION_H |
21 | 21 |
|
22 | 22 |
#include <lemon/tolerance.h> |
23 | 23 |
#include <lemon/elevator.h> |
24 | 24 |
|
25 | 25 |
///\ingroup max_flow |
26 | 26 |
///\file |
27 | 27 |
///\brief Push-relabel algorithm for finding a feasible circulation. |
28 | 28 |
/// |
29 | 29 |
namespace lemon { |
30 | 30 |
|
31 | 31 |
/// \brief Default traits class of Circulation class. |
32 | 32 |
/// |
33 | 33 |
/// Default traits class of Circulation class. |
34 |
/// \tparam GR Digraph type. |
|
35 |
/// \tparam LM Lower bound capacity map type. |
|
36 |
/// \tparam UM Upper bound capacity map type. |
|
37 |
/// \tparam DM Delta map type. |
|
34 |
/// |
|
35 |
/// \tparam GR Type of the digraph the algorithm runs on. |
|
36 |
/// \tparam LM The type of the lower bound map. |
|
37 |
/// \tparam UM The type of the upper bound (capacity) map. |
|
38 |
/// \tparam SM The type of the supply map. |
|
38 | 39 |
template <typename GR, typename LM, |
39 |
typename UM, typename |
|
40 |
typename UM, typename SM> |
|
40 | 41 |
struct CirculationDefaultTraits { |
41 | 42 |
|
42 | 43 |
/// \brief The type of the digraph the algorithm runs on. |
43 | 44 |
typedef GR Digraph; |
44 | 45 |
|
45 |
/// \brief The type of the map that stores the circulation lower |
|
46 |
/// bound. |
|
46 |
/// \brief The type of the lower bound map. |
|
47 | 47 |
/// |
48 |
/// The type of the map that stores the circulation lower bound. |
|
49 |
/// It must meet the \ref concepts::ReadMap "ReadMap" concept. |
|
50 |
|
|
48 |
/// The type of the map that stores the lower bounds on the arcs. |
|
49 |
/// It must conform to the \ref concepts::ReadMap "ReadMap" concept. |
|
50 |
typedef LM LowerMap; |
|
51 | 51 |
|
52 |
/// \brief The type of the map that stores the circulation upper |
|
53 |
/// bound. |
|
52 |
/// \brief The type of the upper bound (capacity) map. |
|
54 | 53 |
/// |
55 |
/// The type of the map that stores the circulation upper bound. |
|
56 |
/// It must meet the \ref concepts::ReadMap "ReadMap" concept. |
|
57 |
|
|
54 |
/// The type of the map that stores the upper bounds (capacities) |
|
55 |
/// on the arcs. |
|
56 |
/// It must conform to the \ref concepts::ReadMap "ReadMap" concept. |
|
57 |
typedef UM UpperMap; |
|
58 | 58 |
|
59 |
/// \brief The type of the map that stores the lower bound for |
|
60 |
/// the supply of the nodes. |
|
59 |
/// \brief The type of supply map. |
|
61 | 60 |
/// |
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 DM DeltaMap; |
|
61 |
/// The type of the map that stores the signed supply values of the |
|
62 |
/// nodes. |
|
63 |
/// It must conform to the \ref concepts::ReadMap "ReadMap" concept. |
|
64 |
typedef SM SupplyMap; |
|
66 | 65 |
|
67 | 66 |
/// \brief The type of the flow values. |
68 |
typedef typename |
|
67 |
typedef typename SupplyMap::Value Flow; |
|
69 | 68 |
|
70 | 69 |
/// \brief The type of the map that stores the flow values. |
71 | 70 |
/// |
72 | 71 |
/// 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; |
|
72 |
/// It must conform to the \ref concepts::ReadWriteMap "ReadWriteMap" |
|
73 |
/// concept. |
|
74 |
typedef typename Digraph::template ArcMap<Flow> FlowMap; |
|
75 | 75 |
|
76 | 76 |
/// \brief Instantiates a FlowMap. |
77 | 77 |
/// |
78 | 78 |
/// This function instantiates a \ref FlowMap. |
79 |
/// \param digraph The digraph |
|
79 |
/// \param digraph The digraph for which we would like to define |
|
80 | 80 |
/// the flow map. |
81 | 81 |
static FlowMap* createFlowMap(const Digraph& digraph) { |
82 | 82 |
return new FlowMap(digraph); |
83 | 83 |
} |
84 | 84 |
|
85 | 85 |
/// \brief The elevator type used by the algorithm. |
86 | 86 |
/// |
87 | 87 |
/// The elevator type used by the algorithm. |
88 | 88 |
/// |
89 | 89 |
/// \sa Elevator |
90 | 90 |
/// \sa LinkedElevator |
91 | 91 |
typedef lemon::Elevator<Digraph, typename Digraph::Node> Elevator; |
92 | 92 |
|
93 | 93 |
/// \brief Instantiates an Elevator. |
94 | 94 |
/// |
95 | 95 |
/// This function instantiates an \ref Elevator. |
96 |
/// \param digraph The digraph |
|
96 |
/// \param digraph The digraph for which we would like to define |
|
97 | 97 |
/// the elevator. |
98 | 98 |
/// \param max_level The maximum level of the elevator. |
99 | 99 |
static Elevator* createElevator(const Digraph& digraph, int max_level) { |
100 | 100 |
return new Elevator(digraph, max_level); |
101 | 101 |
} |
102 | 102 |
|
103 | 103 |
/// \brief The tolerance used by the algorithm |
104 | 104 |
/// |
105 | 105 |
/// The tolerance used by the algorithm to handle inexact computation. |
106 |
typedef lemon::Tolerance< |
|
106 |
typedef lemon::Tolerance<Flow> Tolerance; |
|
107 | 107 |
|
108 | 108 |
}; |
109 | 109 |
|
110 | 110 |
/** |
111 | 111 |
\brief Push-relabel algorithm for the network circulation problem. |
112 | 112 |
|
113 | 113 |
\ingroup max_flow |
114 |
This class implements a push-relabel algorithm for the network |
|
115 |
circulation problem. |
|
114 |
This class implements a push-relabel algorithm for the \e network |
|
115 |
\e circulation problem. |
|
116 | 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. |
|
117 |
are given for the flow values on the arcs and lower bounds are |
|
118 |
given for the difference between the outgoing and incoming flow |
|
119 |
at the nodes. |
|
119 | 120 |
|
120 | 121 |
The exact formulation of this problem is the following. |
121 | 122 |
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. |
|
123 |
\f$lower, upper: A\rightarrow\mathbf{R}^+_0\f$ denote the lower and |
|
124 |
upper bounds on the arcs, for which \f$0 \leq lower(uv) \leq upper(uv)\f$ |
|
125 |
holds for all \f$uv\in A\f$, and \f$sup: V\rightarrow\mathbf{R}\f$ |
|
126 |
denotes the signed supply values of the nodes. |
|
127 |
If \f$sup(u)>0\f$, then \f$u\f$ is a supply node with \f$sup(u)\f$ |
|
128 |
supply, if \f$sup(u)<0\f$, then \f$u\f$ is a demand node with |
|
129 |
\f$-sup(u)\f$ demand. |
|
130 |
A feasible circulation is an \f$f: A\rightarrow\mathbf{R}^+_0\f$ |
|
131 |
solution of the following problem. |
|
132 | 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. |
|
133 |
\f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) |
|
134 |
\geq sup(u) \quad \forall u\in V, \f] |
|
135 |
\f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A. \f] |
|
136 |
|
|
137 |
The sum of the supply values, i.e. \f$\sum_{u\in V} sup(u)\f$ must be |
|
138 |
zero or negative in order to have a feasible solution (since the sum |
|
139 |
of the expressions on the left-hand side of the inequalities is zero). |
|
140 |
It means that the total demand must be greater or equal to the total |
|
141 |
supply and all the supplies have to be carried out from the supply nodes, |
|
142 |
but there could be demands that are not satisfied. |
|
143 |
If \f$\sum_{u\in V} sup(u)\f$ is zero, then all the supply/demand |
|
144 |
constraints have to be satisfied with equality, i.e. all demands |
|
145 |
have to be satisfied and all supplies have to be used. |
|
146 |
|
|
147 |
If you need the opposite inequalities in the supply/demand constraints |
|
148 |
(i.e. the total demand is less than the total supply and all the demands |
|
149 |
have to be satisfied while there could be supplies that are not used), |
|
150 |
then you could easily transform the problem to the above form by reversing |
|
151 |
the direction of the arcs and taking the negative of the supply values |
|
152 |
(e.g. using \ref ReverseDigraph and \ref NegMap adaptors). |
|
153 |
|
|
154 |
Note that this algorithm also provides a feasible solution for the |
|
155 |
\ref min_cost_flow "minimum cost flow problem". |
|
139 | 156 |
|
140 | 157 |
\tparam GR The type of the digraph the algorithm runs on. |
141 |
\tparam LM The type of the lower bound |
|
158 |
\tparam LM The type of the lower bound map. The default |
|
142 | 159 |
map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>". |
143 |
\tparam UM The type of the upper bound capacity map. The default |
|
144 |
map type is \c LM. |
|
145 |
\tparam DM The type of the map that stores the lower bound |
|
146 |
for the supply of the nodes. The default map type is |
|
160 |
\tparam UM The type of the upper bound (capacity) map. |
|
161 |
The default map type is \c LM. |
|
162 |
\tparam SM The type of the supply map. The default map type is |
|
147 | 163 |
\ref concepts::Digraph::NodeMap "GR::NodeMap<UM::Value>". |
148 | 164 |
*/ |
149 | 165 |
#ifdef DOXYGEN |
150 | 166 |
template< typename GR, |
151 | 167 |
typename LM, |
152 | 168 |
typename UM, |
153 |
typename |
|
169 |
typename SM, |
|
154 | 170 |
typename TR > |
155 | 171 |
#else |
156 | 172 |
template< typename GR, |
157 | 173 |
typename LM = typename GR::template ArcMap<int>, |
158 | 174 |
typename UM = LM, |
159 |
typename DM = typename GR::template NodeMap<typename UM::Value>, |
|
160 |
typename TR = CirculationDefaultTraits<GR, LM, UM, DM> > |
|
175 |
typename SM = typename GR::template NodeMap<typename UM::Value>, |
|
176 |
typename TR = CirculationDefaultTraits<GR, LM, UM, SM> > |
|
161 | 177 |
#endif |
162 | 178 |
class Circulation { |
163 | 179 |
public: |
164 | 180 |
|
165 | 181 |
///The \ref CirculationDefaultTraits "traits class" of the algorithm. |
166 | 182 |
typedef TR Traits; |
167 | 183 |
///The type of the digraph the algorithm runs on. |
168 | 184 |
typedef typename Traits::Digraph Digraph; |
169 | 185 |
///The type of the flow values. |
170 |
typedef typename Traits:: |
|
186 |
typedef typename Traits::Flow Flow; |
|
171 | 187 |
|
172 |
/// The type of the lower bound capacity map. |
|
173 |
typedef typename Traits::LCapMap LCapMap; |
|
174 |
/// The type of the upper bound capacity map. |
|
175 |
typedef typename Traits::UCapMap UCapMap; |
|
176 |
/// \brief The type of the map that stores the lower bound for |
|
177 |
/// the supply of the nodes. |
|
178 |
|
|
188 |
///The type of the lower bound map. |
|
189 |
typedef typename Traits::LowerMap LowerMap; |
|
190 |
///The type of the upper bound (capacity) map. |
|
191 |
typedef typename Traits::UpperMap UpperMap; |
|
192 |
///The type of the supply map. |
|
193 |
typedef typename Traits::SupplyMap SupplyMap; |
|
179 | 194 |
///The type of the flow map. |
180 | 195 |
typedef typename Traits::FlowMap FlowMap; |
181 | 196 |
|
182 | 197 |
///The type of the elevator. |
183 | 198 |
typedef typename Traits::Elevator Elevator; |
184 | 199 |
///The type of the tolerance. |
185 | 200 |
typedef typename Traits::Tolerance Tolerance; |
186 | 201 |
|
187 | 202 |
private: |
188 | 203 |
|
189 | 204 |
TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); |
190 | 205 |
|
191 | 206 |
const Digraph &_g; |
192 | 207 |
int _node_num; |
193 | 208 |
|
194 |
const LCapMap *_lo; |
|
195 |
const UCapMap *_up; |
|
196 |
const |
|
209 |
const LowerMap *_lo; |
|
210 |
const UpperMap *_up; |
|
211 |
const SupplyMap *_supply; |
|
197 | 212 |
|
198 | 213 |
FlowMap *_flow; |
199 | 214 |
bool _local_flow; |
200 | 215 |
|
201 | 216 |
Elevator* _level; |
202 | 217 |
bool _local_level; |
203 | 218 |
|
204 |
typedef typename Digraph::template NodeMap< |
|
219 |
typedef typename Digraph::template NodeMap<Flow> ExcessMap; |
|
205 | 220 |
ExcessMap* _excess; |
206 | 221 |
|
207 | 222 |
Tolerance _tol; |
208 | 223 |
int _el; |
209 | 224 |
|
210 | 225 |
public: |
211 | 226 |
|
212 | 227 |
typedef Circulation Create; |
213 | 228 |
|
214 | 229 |
///\name Named Template Parameters |
215 | 230 |
|
216 | 231 |
///@{ |
217 | 232 |
|
218 | 233 |
template <typename _FlowMap> |
219 | 234 |
struct SetFlowMapTraits : public Traits { |
220 | 235 |
typedef _FlowMap FlowMap; |
221 | 236 |
static FlowMap *createFlowMap(const Digraph&) { |
222 | 237 |
LEMON_ASSERT(false, "FlowMap is not initialized"); |
223 | 238 |
return 0; // ignore warnings |
224 | 239 |
} |
225 | 240 |
}; |
226 | 241 |
|
227 | 242 |
/// \brief \ref named-templ-param "Named parameter" for setting |
228 | 243 |
/// FlowMap type |
229 | 244 |
/// |
230 | 245 |
/// \ref named-templ-param "Named parameter" for setting FlowMap |
231 | 246 |
/// type. |
232 | 247 |
template <typename _FlowMap> |
233 | 248 |
struct SetFlowMap |
234 |
: public Circulation<Digraph, |
|
249 |
: public Circulation<Digraph, LowerMap, UpperMap, SupplyMap, |
|
235 | 250 |
SetFlowMapTraits<_FlowMap> > { |
236 |
typedef Circulation<Digraph, |
|
251 |
typedef Circulation<Digraph, LowerMap, UpperMap, SupplyMap, |
|
237 | 252 |
SetFlowMapTraits<_FlowMap> > Create; |
238 | 253 |
}; |
239 | 254 |
|
240 | 255 |
template <typename _Elevator> |
241 | 256 |
struct SetElevatorTraits : public Traits { |
242 | 257 |
typedef _Elevator Elevator; |
243 | 258 |
static Elevator *createElevator(const Digraph&, int) { |
244 | 259 |
LEMON_ASSERT(false, "Elevator is not initialized"); |
245 | 260 |
return 0; // ignore warnings |
246 | 261 |
} |
247 | 262 |
}; |
248 | 263 |
|
249 | 264 |
/// \brief \ref named-templ-param "Named parameter" for setting |
250 | 265 |
/// Elevator type |
251 | 266 |
/// |
252 | 267 |
/// \ref named-templ-param "Named parameter" for setting Elevator |
253 | 268 |
/// type. If this named parameter is used, then an external |
254 | 269 |
/// elevator object must be passed to the algorithm using the |
255 | 270 |
/// \ref elevator(Elevator&) "elevator()" function before calling |
256 | 271 |
/// \ref run() or \ref init(). |
257 | 272 |
/// \sa SetStandardElevator |
258 | 273 |
template <typename _Elevator> |
259 | 274 |
struct SetElevator |
260 |
: public Circulation<Digraph, |
|
275 |
: public Circulation<Digraph, LowerMap, UpperMap, SupplyMap, |
|
261 | 276 |
SetElevatorTraits<_Elevator> > { |
262 |
typedef Circulation<Digraph, |
|
277 |
typedef Circulation<Digraph, LowerMap, UpperMap, SupplyMap, |
|
263 | 278 |
SetElevatorTraits<_Elevator> > Create; |
264 | 279 |
}; |
265 | 280 |
|
266 | 281 |
template <typename _Elevator> |
267 | 282 |
struct SetStandardElevatorTraits : public Traits { |
268 | 283 |
typedef _Elevator Elevator; |
269 | 284 |
static Elevator *createElevator(const Digraph& digraph, int max_level) { |
270 | 285 |
return new Elevator(digraph, max_level); |
271 | 286 |
} |
272 | 287 |
}; |
273 | 288 |
|
274 | 289 |
/// \brief \ref named-templ-param "Named parameter" for setting |
275 | 290 |
/// Elevator type with automatic allocation |
276 | 291 |
/// |
277 | 292 |
/// \ref named-templ-param "Named parameter" for setting Elevator |
278 | 293 |
/// type with automatic allocation. |
279 | 294 |
/// The Elevator should have standard constructor interface to be |
280 | 295 |
/// able to automatically created by the algorithm (i.e. the |
281 | 296 |
/// digraph and the maximum level should be passed to it). |
282 | 297 |
/// However an external elevator object could also be passed to the |
283 | 298 |
/// algorithm with the \ref elevator(Elevator&) "elevator()" function |
284 | 299 |
/// before calling \ref run() or \ref init(). |
285 | 300 |
/// \sa SetElevator |
286 | 301 |
template <typename _Elevator> |
287 | 302 |
struct SetStandardElevator |
288 |
: public Circulation<Digraph, |
|
303 |
: public Circulation<Digraph, LowerMap, UpperMap, SupplyMap, |
|
289 | 304 |
SetStandardElevatorTraits<_Elevator> > { |
290 |
typedef Circulation<Digraph, |
|
305 |
typedef Circulation<Digraph, LowerMap, UpperMap, SupplyMap, |
|
291 | 306 |
SetStandardElevatorTraits<_Elevator> > Create; |
292 | 307 |
}; |
293 | 308 |
|
294 | 309 |
/// @} |
295 | 310 |
|
296 | 311 |
protected: |
297 | 312 |
|
298 | 313 |
Circulation() {} |
299 | 314 |
|
300 | 315 |
public: |
301 | 316 |
|
302 |
/// |
|
317 |
/// Constructor. |
|
303 | 318 |
|
304 | 319 |
/// The constructor of the class. |
305 |
/// \param g The digraph the algorithm runs on. |
|
306 |
/// \param lo The lower bound capacity of the arcs. |
|
307 |
/// \param up The upper bound capacity of the arcs. |
|
308 |
/// \param delta The lower bound for the supply of the nodes. |
|
309 |
Circulation(const Digraph &g,const LCapMap &lo, |
|
310 |
const UCapMap &up,const DeltaMap &delta) |
|
311 |
: _g(g), _node_num(), |
|
312 |
_lo(&lo),_up(&up),_delta(&delta),_flow(0),_local_flow(false), |
|
313 |
|
|
320 |
/// |
|
321 |
/// \param graph The digraph the algorithm runs on. |
|
322 |
/// \param lower The lower bounds for the flow values on the arcs. |
|
323 |
/// \param upper The upper bounds (capacities) for the flow values |
|
324 |
/// on the arcs. |
|
325 |
/// \param supply The signed supply values of the nodes. |
|
326 |
Circulation(const Digraph &graph, const LowerMap &lower, |
|
327 |
const UpperMap &upper, const SupplyMap &supply) |
|
328 |
: _g(graph), _lo(&lower), _up(&upper), _supply(&supply), |
|
329 |
_flow(NULL), _local_flow(false), _level(NULL), _local_level(false), |
|
330 |
_excess(NULL) {} |
|
314 | 331 |
|
315 | 332 |
/// Destructor. |
316 | 333 |
~Circulation() { |
317 | 334 |
destroyStructures(); |
318 | 335 |
} |
319 | 336 |
|
320 | 337 |
|
321 | 338 |
private: |
322 | 339 |
|
323 | 340 |
void createStructures() { |
324 | 341 |
_node_num = _el = countNodes(_g); |
325 | 342 |
|
326 | 343 |
if (!_flow) { |
327 | 344 |
_flow = Traits::createFlowMap(_g); |
328 | 345 |
_local_flow = true; |
329 | 346 |
} |
330 | 347 |
if (!_level) { |
331 | 348 |
_level = Traits::createElevator(_g, _node_num); |
332 | 349 |
_local_level = true; |
333 | 350 |
} |
334 | 351 |
if (!_excess) { |
335 | 352 |
_excess = new ExcessMap(_g); |
336 | 353 |
} |
337 | 354 |
} |
338 | 355 |
|
339 | 356 |
void destroyStructures() { |
340 | 357 |
if (_local_flow) { |
341 | 358 |
delete _flow; |
342 | 359 |
} |
343 | 360 |
if (_local_level) { |
344 | 361 |
delete _level; |
345 | 362 |
} |
346 | 363 |
if (_excess) { |
347 | 364 |
delete _excess; |
348 | 365 |
} |
349 | 366 |
} |
350 | 367 |
|
351 | 368 |
public: |
352 | 369 |
|
353 |
/// Sets the lower bound |
|
370 |
/// Sets the lower bound map. |
|
354 | 371 |
|
355 |
/// Sets the lower bound |
|
372 |
/// Sets the lower bound map. |
|
356 | 373 |
/// \return <tt>(*this)</tt> |
357 |
Circulation& |
|
374 |
Circulation& lowerMap(const LowerMap& map) { |
|
358 | 375 |
_lo = ↦ |
359 | 376 |
return *this; |
360 | 377 |
} |
361 | 378 |
|
362 |
/// Sets the upper bound capacity map. |
|
379 |
/// Sets the upper bound (capacity) map. |
|
363 | 380 |
|
364 |
/// Sets the upper bound capacity map. |
|
381 |
/// Sets the upper bound (capacity) map. |
|
365 | 382 |
/// \return <tt>(*this)</tt> |
366 |
Circulation& |
|
383 |
Circulation& upperMap(const LowerMap& map) { |
|
367 | 384 |
_up = ↦ |
368 | 385 |
return *this; |
369 | 386 |
} |
370 | 387 |
|
371 |
/// Sets the |
|
388 |
/// Sets the supply map. |
|
372 | 389 |
|
373 |
/// Sets the |
|
390 |
/// Sets the supply map. |
|
374 | 391 |
/// \return <tt>(*this)</tt> |
375 |
Circulation& deltaMap(const DeltaMap& map) { |
|
376 |
_delta = ↦ |
|
392 |
Circulation& supplyMap(const SupplyMap& map) { |
|
393 |
_supply = ↦ |
|
377 | 394 |
return *this; |
378 | 395 |
} |
379 | 396 |
|
380 | 397 |
/// \brief Sets the flow map. |
381 | 398 |
/// |
382 | 399 |
/// Sets the flow map. |
383 | 400 |
/// If you don't use this function before calling \ref run() or |
384 | 401 |
/// \ref init(), an instance will be allocated automatically. |
385 | 402 |
/// The destructor deallocates this automatically allocated map, |
386 | 403 |
/// of course. |
387 | 404 |
/// \return <tt>(*this)</tt> |
388 | 405 |
Circulation& flowMap(FlowMap& map) { |
389 | 406 |
if (_local_flow) { |
390 | 407 |
delete _flow; |
391 | 408 |
_local_flow = false; |
392 | 409 |
} |
393 | 410 |
_flow = ↦ |
394 | 411 |
return *this; |
395 | 412 |
} |
396 | 413 |
|
397 | 414 |
/// \brief Sets the elevator used by algorithm. |
398 | 415 |
/// |
399 | 416 |
/// Sets the elevator used by algorithm. |
400 | 417 |
/// If you don't use this function before calling \ref run() or |
401 | 418 |
/// \ref init(), an instance will be allocated automatically. |
402 | 419 |
/// The destructor deallocates this automatically allocated elevator, |
403 | 420 |
/// of course. |
404 | 421 |
/// \return <tt>(*this)</tt> |
405 | 422 |
Circulation& elevator(Elevator& elevator) { |
406 | 423 |
if (_local_level) { |
407 | 424 |
delete _level; |
408 | 425 |
_local_level = false; |
409 | 426 |
} |
410 | 427 |
_level = &elevator; |
411 | 428 |
return *this; |
412 | 429 |
} |
413 | 430 |
|
414 | 431 |
/// \brief Returns a const reference to the elevator. |
415 | 432 |
/// |
416 | 433 |
/// Returns a const reference to the elevator. |
417 | 434 |
/// |
418 | 435 |
/// \pre Either \ref run() or \ref init() must be called before |
419 | 436 |
/// using this function. |
420 | 437 |
const Elevator& elevator() const { |
421 | 438 |
return *_level; |
422 | 439 |
} |
423 | 440 |
|
424 | 441 |
/// \brief Sets the tolerance used by algorithm. |
425 | 442 |
/// |
426 | 443 |
/// Sets the tolerance used by algorithm. |
427 | 444 |
Circulation& tolerance(const Tolerance& tolerance) const { |
428 | 445 |
_tol = tolerance; |
429 | 446 |
return *this; |
430 | 447 |
} |
431 | 448 |
|
432 | 449 |
/// \brief Returns a const reference to the tolerance. |
433 | 450 |
/// |
434 | 451 |
/// Returns a const reference to the tolerance. |
435 | 452 |
const Tolerance& tolerance() const { |
436 | 453 |
return tolerance; |
437 | 454 |
} |
438 | 455 |
|
439 | 456 |
/// \name Execution Control |
440 | 457 |
/// The simplest way to execute the algorithm is to call \ref run().\n |
441 | 458 |
/// If you need more control on the initial solution or the execution, |
442 | 459 |
/// first you have to call one of the \ref init() functions, then |
443 | 460 |
/// the \ref start() function. |
444 | 461 |
|
445 | 462 |
///@{ |
446 | 463 |
|
447 | 464 |
/// Initializes the internal data structures. |
448 | 465 |
|
449 | 466 |
/// Initializes the internal data structures and sets all flow values |
450 | 467 |
/// to the lower bound. |
451 | 468 |
void init() |
452 | 469 |
{ |
453 | 470 |
createStructures(); |
454 | 471 |
|
455 | 472 |
for(NodeIt n(_g);n!=INVALID;++n) { |
456 |
_excess->set(n, (* |
|
473 |
_excess->set(n, (*_supply)[n]); |
|
457 | 474 |
} |
458 | 475 |
|
459 | 476 |
for (ArcIt e(_g);e!=INVALID;++e) { |
460 | 477 |
_flow->set(e, (*_lo)[e]); |
461 | 478 |
_excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_flow)[e]); |
462 | 479 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_flow)[e]); |
463 | 480 |
} |
464 | 481 |
|
465 | 482 |
// global relabeling tested, but in general case it provides |
466 | 483 |
// worse performance for random digraphs |
467 | 484 |
_level->initStart(); |
468 | 485 |
for(NodeIt n(_g);n!=INVALID;++n) |
469 | 486 |
_level->initAddItem(n); |
470 | 487 |
_level->initFinish(); |
471 | 488 |
for(NodeIt n(_g);n!=INVALID;++n) |
472 | 489 |
if(_tol.positive((*_excess)[n])) |
473 | 490 |
_level->activate(n); |
474 | 491 |
} |
475 | 492 |
|
476 | 493 |
/// Initializes the internal data structures using a greedy approach. |
477 | 494 |
|
478 | 495 |
/// Initializes the internal data structures using a greedy approach |
479 | 496 |
/// to construct the initial solution. |
480 | 497 |
void greedyInit() |
481 | 498 |
{ |
482 | 499 |
createStructures(); |
483 | 500 |
|
484 | 501 |
for(NodeIt n(_g);n!=INVALID;++n) { |
485 |
_excess->set(n, (* |
|
502 |
_excess->set(n, (*_supply)[n]); |
|
486 | 503 |
} |
487 | 504 |
|
488 | 505 |
for (ArcIt e(_g);e!=INVALID;++e) { |
489 | 506 |
if (!_tol.positive((*_excess)[_g.target(e)] + (*_up)[e])) { |
490 | 507 |
_flow->set(e, (*_up)[e]); |
491 | 508 |
_excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_up)[e]); |
492 | 509 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_up)[e]); |
493 | 510 |
} else if (_tol.positive((*_excess)[_g.target(e)] + (*_lo)[e])) { |
494 | 511 |
_flow->set(e, (*_lo)[e]); |
495 | 512 |
_excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_lo)[e]); |
496 | 513 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_lo)[e]); |
497 | 514 |
} else { |
498 |
|
|
515 |
Flow fc = -(*_excess)[_g.target(e)]; |
|
499 | 516 |
_flow->set(e, fc); |
500 | 517 |
_excess->set(_g.target(e), 0); |
501 | 518 |
_excess->set(_g.source(e), (*_excess)[_g.source(e)] - fc); |
502 | 519 |
} |
503 | 520 |
} |
504 | 521 |
|
505 | 522 |
_level->initStart(); |
506 | 523 |
for(NodeIt n(_g);n!=INVALID;++n) |
507 | 524 |
_level->initAddItem(n); |
508 | 525 |
_level->initFinish(); |
509 | 526 |
for(NodeIt n(_g);n!=INVALID;++n) |
510 | 527 |
if(_tol.positive((*_excess)[n])) |
511 | 528 |
_level->activate(n); |
512 | 529 |
} |
513 | 530 |
|
514 | 531 |
///Executes the algorithm |
515 | 532 |
|
516 | 533 |
///This function executes the algorithm. |
517 | 534 |
/// |
518 | 535 |
///\return \c true if a feasible circulation is found. |
519 | 536 |
/// |
520 | 537 |
///\sa barrier() |
521 | 538 |
///\sa barrierMap() |
522 | 539 |
bool start() |
523 | 540 |
{ |
524 | 541 |
|
525 | 542 |
Node act; |
526 | 543 |
Node bact=INVALID; |
527 | 544 |
Node last_activated=INVALID; |
528 | 545 |
while((act=_level->highestActive())!=INVALID) { |
529 | 546 |
int actlevel=(*_level)[act]; |
530 | 547 |
int mlevel=_node_num; |
531 |
|
|
548 |
Flow exc=(*_excess)[act]; |
|
532 | 549 |
|
533 | 550 |
for(OutArcIt e(_g,act);e!=INVALID; ++e) { |
534 | 551 |
Node v = _g.target(e); |
535 |
|
|
552 |
Flow fc=(*_up)[e]-(*_flow)[e]; |
|
536 | 553 |
if(!_tol.positive(fc)) continue; |
537 | 554 |
if((*_level)[v]<actlevel) { |
538 | 555 |
if(!_tol.less(fc, exc)) { |
539 | 556 |
_flow->set(e, (*_flow)[e] + exc); |
540 | 557 |
_excess->set(v, (*_excess)[v] + exc); |
541 | 558 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
542 | 559 |
_level->activate(v); |
543 | 560 |
_excess->set(act,0); |
544 | 561 |
_level->deactivate(act); |
545 | 562 |
goto next_l; |
546 | 563 |
} |
547 | 564 |
else { |
548 | 565 |
_flow->set(e, (*_up)[e]); |
549 | 566 |
_excess->set(v, (*_excess)[v] + fc); |
550 | 567 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
551 | 568 |
_level->activate(v); |
552 | 569 |
exc-=fc; |
553 | 570 |
} |
554 | 571 |
} |
555 | 572 |
else if((*_level)[v]<mlevel) mlevel=(*_level)[v]; |
556 | 573 |
} |
557 | 574 |
for(InArcIt e(_g,act);e!=INVALID; ++e) { |
558 | 575 |
Node v = _g.source(e); |
559 |
|
|
576 |
Flow fc=(*_flow)[e]-(*_lo)[e]; |
|
560 | 577 |
if(!_tol.positive(fc)) continue; |
561 | 578 |
if((*_level)[v]<actlevel) { |
562 | 579 |
if(!_tol.less(fc, exc)) { |
563 | 580 |
_flow->set(e, (*_flow)[e] - exc); |
564 | 581 |
_excess->set(v, (*_excess)[v] + exc); |
565 | 582 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
566 | 583 |
_level->activate(v); |
567 | 584 |
_excess->set(act,0); |
568 | 585 |
_level->deactivate(act); |
569 | 586 |
goto next_l; |
570 | 587 |
} |
571 | 588 |
else { |
572 | 589 |
_flow->set(e, (*_lo)[e]); |
573 | 590 |
_excess->set(v, (*_excess)[v] + fc); |
574 | 591 |
if(!_level->active(v) && _tol.positive((*_excess)[v])) |
575 | 592 |
_level->activate(v); |
576 | 593 |
exc-=fc; |
577 | 594 |
} |
578 | 595 |
} |
579 | 596 |
else if((*_level)[v]<mlevel) mlevel=(*_level)[v]; |
580 | 597 |
} |
581 | 598 |
|
582 | 599 |
_excess->set(act, exc); |
583 | 600 |
if(!_tol.positive(exc)) _level->deactivate(act); |
584 | 601 |
else if(mlevel==_node_num) { |
585 | 602 |
_level->liftHighestActiveToTop(); |
586 | 603 |
_el = _node_num; |
587 | 604 |
return false; |
588 | 605 |
} |
589 | 606 |
else { |
590 | 607 |
_level->liftHighestActive(mlevel+1); |
591 | 608 |
if(_level->onLevel(actlevel)==0) { |
592 | 609 |
_el = actlevel; |
593 | 610 |
return false; |
594 | 611 |
} |
595 | 612 |
} |
596 | 613 |
next_l: |
597 | 614 |
; |
598 | 615 |
} |
599 | 616 |
return true; |
600 | 617 |
} |
601 | 618 |
|
602 | 619 |
/// Runs the algorithm. |
603 | 620 |
|
604 | 621 |
/// This function runs the algorithm. |
605 | 622 |
/// |
606 | 623 |
/// \return \c true if a feasible circulation is found. |
607 | 624 |
/// |
608 | 625 |
/// \note Apart from the return value, c.run() is just a shortcut of |
609 | 626 |
/// the following code. |
610 | 627 |
/// \code |
611 | 628 |
/// c.greedyInit(); |
612 | 629 |
/// c.start(); |
613 | 630 |
/// \endcode |
614 | 631 |
bool run() { |
615 | 632 |
greedyInit(); |
616 | 633 |
return start(); |
617 | 634 |
} |
618 | 635 |
|
619 | 636 |
/// @} |
620 | 637 |
|
621 | 638 |
/// \name Query Functions |
622 | 639 |
/// The results of the circulation algorithm can be obtained using |
623 | 640 |
/// these functions.\n |
624 | 641 |
/// Either \ref run() or \ref start() should be called before |
625 | 642 |
/// using them. |
626 | 643 |
|
627 | 644 |
///@{ |
628 | 645 |
|
629 | 646 |
/// \brief Returns the flow on the given arc. |
630 | 647 |
/// |
631 | 648 |
/// Returns the flow on the given arc. |
632 | 649 |
/// |
633 | 650 |
/// \pre Either \ref run() or \ref init() must be called before |
634 | 651 |
/// using this function. |
635 |
|
|
652 |
Flow flow(const Arc& arc) const { |
|
636 | 653 |
return (*_flow)[arc]; |
637 | 654 |
} |
638 | 655 |
|
639 | 656 |
/// \brief Returns a const reference to the flow map. |
640 | 657 |
/// |
641 | 658 |
/// Returns a const reference to the arc map storing the found flow. |
642 | 659 |
/// |
643 | 660 |
/// \pre Either \ref run() or \ref init() must be called before |
644 | 661 |
/// using this function. |
645 | 662 |
const FlowMap& flowMap() const { |
646 | 663 |
return *_flow; |
647 | 664 |
} |
648 | 665 |
|
649 | 666 |
/** |
650 | 667 |
\brief Returns \c true if the given node is in a barrier. |
651 | 668 |
|
652 | 669 |
Barrier is a set \e B of nodes for which |
653 | 670 |
|
654 |
\f[ \sum_{a\in\delta_{out}(B)} upper(a) - |
|
655 |
\sum_{a\in\delta_{in}(B)} lower(a) < \sum_{v\in B}delta(v) \f] |
|
671 |
\f[ \sum_{uv\in A: u\in B} upper(uv) - |
|
672 |
\sum_{uv\in A: v\in B} lower(uv) < \sum_{v\in B} sup(v) \f] |
|
656 | 673 |
|
657 | 674 |
holds. The existence of a set with this property prooves that a |
658 | 675 |
feasible circualtion cannot exist. |
659 | 676 |
|
660 | 677 |
This function returns \c true if the given node is in the found |
661 | 678 |
barrier. If a feasible circulation is found, the function |
662 | 679 |
gives back \c false for every node. |
663 | 680 |
|
664 | 681 |
\pre Either \ref run() or \ref init() must be called before |
665 | 682 |
using this function. |
666 | 683 |
|
667 | 684 |
\sa barrierMap() |
668 | 685 |
\sa checkBarrier() |
669 | 686 |
*/ |
670 | 687 |
bool barrier(const Node& node) const |
671 | 688 |
{ |
672 | 689 |
return (*_level)[node] >= _el; |
673 | 690 |
} |
674 | 691 |
|
675 | 692 |
/// \brief Gives back a barrier. |
676 | 693 |
/// |
677 | 694 |
/// This function sets \c bar to the characteristic vector of the |
678 | 695 |
/// found barrier. \c bar should be a \ref concepts::WriteMap "writable" |
679 | 696 |
/// node map with \c bool (or convertible) value type. |
680 | 697 |
/// |
681 | 698 |
/// If a feasible circulation is found, the function gives back an |
682 | 699 |
/// empty set, so \c bar[v] will be \c false for all nodes \c v. |
683 | 700 |
/// |
684 | 701 |
/// \note This function calls \ref barrier() for each node, |
685 | 702 |
/// so it runs in \f$O(n)\f$ time. |
686 | 703 |
/// |
687 | 704 |
/// \pre Either \ref run() or \ref init() must be called before |
688 | 705 |
/// using this function. |
689 | 706 |
/// |
690 | 707 |
/// \sa barrier() |
691 | 708 |
/// \sa checkBarrier() |
692 | 709 |
template<class BarrierMap> |
693 | 710 |
void barrierMap(BarrierMap &bar) const |
694 | 711 |
{ |
695 | 712 |
for(NodeIt n(_g);n!=INVALID;++n) |
696 | 713 |
bar.set(n, (*_level)[n] >= _el); |
697 | 714 |
} |
698 | 715 |
|
699 | 716 |
/// @} |
700 | 717 |
|
701 | 718 |
/// \name Checker Functions |
702 | 719 |
/// The feasibility of the results can be checked using |
703 | 720 |
/// these functions.\n |
704 | 721 |
/// Either \ref run() or \ref start() should be called before |
705 | 722 |
/// using them. |
706 | 723 |
|
707 | 724 |
///@{ |
708 | 725 |
|
709 | 726 |
///Check if the found flow is a feasible circulation |
710 | 727 |
|
711 | 728 |
///Check if the found flow is a feasible circulation, |
712 | 729 |
/// |
713 | 730 |
bool checkFlow() const { |
714 | 731 |
for(ArcIt e(_g);e!=INVALID;++e) |
715 | 732 |
if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false; |
716 | 733 |
for(NodeIt n(_g);n!=INVALID;++n) |
717 | 734 |
{ |
718 |
|
|
735 |
Flow dif=-(*_supply)[n]; |
|
719 | 736 |
for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e]; |
720 | 737 |
for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e]; |
721 | 738 |
if(_tol.negative(dif)) return false; |
722 | 739 |
} |
723 | 740 |
return true; |
724 | 741 |
} |
725 | 742 |
|
726 | 743 |
///Check whether or not the last execution provides a barrier |
727 | 744 |
|
728 | 745 |
///Check whether or not the last execution provides a barrier. |
729 | 746 |
///\sa barrier() |
730 | 747 |
///\sa barrierMap() |
731 | 748 |
bool checkBarrier() const |
732 | 749 |
{ |
733 |
|
|
750 |
Flow delta=0; |
|
734 | 751 |
for(NodeIt n(_g);n!=INVALID;++n) |
735 | 752 |
if(barrier(n)) |
736 |
delta-=(* |
|
753 |
delta-=(*_supply)[n]; |
|
737 | 754 |
for(ArcIt e(_g);e!=INVALID;++e) |
738 | 755 |
{ |
739 | 756 |
Node s=_g.source(e); |
740 | 757 |
Node t=_g.target(e); |
741 | 758 |
if(barrier(s)&&!barrier(t)) delta+=(*_up)[e]; |
742 | 759 |
else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e]; |
743 | 760 |
} |
744 | 761 |
return _tol.negative(delta); |
745 | 762 |
} |
746 | 763 |
|
747 | 764 |
/// @} |
748 | 765 |
|
749 | 766 |
}; |
750 | 767 |
|
751 | 768 |
} |
752 | 769 |
|
753 | 770 |
#endif |
1 | 1 |
/* -*- mode: C++; indent-tabs-mode: nil; -*- |
2 | 2 |
* |
3 | 3 |
* This file is a part of LEMON, a generic C++ optimization library. |
4 | 4 |
* |
5 | 5 |
* Copyright (C) 2003-2009 |
6 | 6 |
* Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport |
7 | 7 |
* (Egervary Research Group on Combinatorial Optimization, EGRES). |
8 | 8 |
* |
9 | 9 |
* Permission to use, modify and distribute this software is granted |
10 | 10 |
* provided that this copyright notice appears in all copies. For |
11 | 11 |
* precise terms see the accompanying LICENSE file. |
12 | 12 |
* |
13 | 13 |
* This software is provided "AS IS" with no warranty of any kind, |
14 | 14 |
* express or implied, and with no claim as to its suitability for any |
15 | 15 |
* purpose. |
16 | 16 |
* |
17 | 17 |
*/ |
18 | 18 |
|
19 | 19 |
#ifndef LEMON_PREFLOW_H |
20 | 20 |
#define LEMON_PREFLOW_H |
21 | 21 |
|
22 | 22 |
#include <lemon/tolerance.h> |
23 | 23 |
#include <lemon/elevator.h> |
24 | 24 |
|
25 | 25 |
/// \file |
26 | 26 |
/// \ingroup max_flow |
27 | 27 |
/// \brief Implementation of the preflow algorithm. |
28 | 28 |
|
29 | 29 |
namespace lemon { |
30 | 30 |
|
31 | 31 |
/// \brief Default traits class of Preflow class. |
32 | 32 |
/// |
33 | 33 |
/// Default traits class of Preflow class. |
34 | 34 |
/// \tparam GR Digraph type. |
35 | 35 |
/// \tparam CM Capacity map type. |
36 | 36 |
template <typename GR, typename CM> |
37 | 37 |
struct PreflowDefaultTraits { |
38 | 38 |
|
39 | 39 |
/// \brief The type of the digraph the algorithm runs on. |
40 | 40 |
typedef GR Digraph; |
41 | 41 |
|
42 | 42 |
/// \brief The type of the map that stores the arc capacities. |
43 | 43 |
/// |
44 | 44 |
/// The type of the map that stores the arc capacities. |
45 | 45 |
/// It must meet the \ref concepts::ReadMap "ReadMap" concept. |
46 | 46 |
typedef CM CapacityMap; |
47 | 47 |
|
48 | 48 |
/// \brief The type of the flow values. |
49 |
typedef typename CapacityMap::Value |
|
49 |
typedef typename CapacityMap::Value Flow; |
|
50 | 50 |
|
51 | 51 |
/// \brief The type of the map that stores the flow values. |
52 | 52 |
/// |
53 | 53 |
/// The type of the map that stores the flow values. |
54 | 54 |
/// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept. |
55 |
typedef typename Digraph::template ArcMap< |
|
55 |
typedef typename Digraph::template ArcMap<Flow> FlowMap; |
|
56 | 56 |
|
57 | 57 |
/// \brief Instantiates a FlowMap. |
58 | 58 |
/// |
59 | 59 |
/// This function instantiates a \ref FlowMap. |
60 |
/// \param digraph The digraph |
|
60 |
/// \param digraph The digraph for which we would like to define |
|
61 | 61 |
/// the flow map. |
62 | 62 |
static FlowMap* createFlowMap(const Digraph& digraph) { |
63 | 63 |
return new FlowMap(digraph); |
64 | 64 |
} |
65 | 65 |
|
66 | 66 |
/// \brief The elevator type used by Preflow algorithm. |
67 | 67 |
/// |
68 | 68 |
/// The elevator type used by Preflow algorithm. |
69 | 69 |
/// |
70 | 70 |
/// \sa Elevator |
71 | 71 |
/// \sa LinkedElevator |
72 | 72 |
typedef LinkedElevator<Digraph, typename Digraph::Node> Elevator; |
73 | 73 |
|
74 | 74 |
/// \brief Instantiates an Elevator. |
75 | 75 |
/// |
76 | 76 |
/// This function instantiates an \ref Elevator. |
77 |
/// \param digraph The digraph |
|
77 |
/// \param digraph The digraph for which we would like to define |
|
78 | 78 |
/// the elevator. |
79 | 79 |
/// \param max_level The maximum level of the elevator. |
80 | 80 |
static Elevator* createElevator(const Digraph& digraph, int max_level) { |
81 | 81 |
return new Elevator(digraph, max_level); |
82 | 82 |
} |
83 | 83 |
|
84 | 84 |
/// \brief The tolerance used by the algorithm |
85 | 85 |
/// |
86 | 86 |
/// The tolerance used by the algorithm to handle inexact computation. |
87 |
typedef lemon::Tolerance< |
|
87 |
typedef lemon::Tolerance<Flow> Tolerance; |
|
88 | 88 |
|
89 | 89 |
}; |
90 | 90 |
|
91 | 91 |
|
92 | 92 |
/// \ingroup max_flow |
93 | 93 |
/// |
94 | 94 |
/// \brief %Preflow algorithm class. |
95 | 95 |
/// |
96 | 96 |
/// This class provides an implementation of Goldberg-Tarjan's \e preflow |
97 | 97 |
/// \e push-relabel algorithm producing a flow of maximum value in a |
98 | 98 |
/// digraph. The preflow algorithms are the fastest known maximum |
99 | 99 |
/// flow algorithms. The current implementation use a mixture of the |
100 | 100 |
/// \e "highest label" and the \e "bound decrease" heuristics. |
101 | 101 |
/// The worst case time complexity of the algorithm is \f$O(n^2\sqrt{e})\f$. |
102 | 102 |
/// |
103 | 103 |
/// The algorithm consists of two phases. After the first phase |
104 | 104 |
/// the maximum flow value and the minimum cut is obtained. The |
105 | 105 |
/// second phase constructs a feasible maximum flow on each arc. |
106 | 106 |
/// |
107 | 107 |
/// \tparam GR The type of the digraph the algorithm runs on. |
108 | 108 |
/// \tparam CM The type of the capacity map. The default map |
109 | 109 |
/// type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>". |
110 | 110 |
#ifdef DOXYGEN |
111 | 111 |
template <typename GR, typename CM, typename TR> |
112 | 112 |
#else |
113 | 113 |
template <typename GR, |
114 | 114 |
typename CM = typename GR::template ArcMap<int>, |
115 | 115 |
typename TR = PreflowDefaultTraits<GR, CM> > |
116 | 116 |
#endif |
117 | 117 |
class Preflow { |
118 | 118 |
public: |
119 | 119 |
|
120 | 120 |
///The \ref PreflowDefaultTraits "traits class" of the algorithm. |
121 | 121 |
typedef TR Traits; |
122 | 122 |
///The type of the digraph the algorithm runs on. |
123 | 123 |
typedef typename Traits::Digraph Digraph; |
124 | 124 |
///The type of the capacity map. |
125 | 125 |
typedef typename Traits::CapacityMap CapacityMap; |
126 | 126 |
///The type of the flow values. |
127 |
typedef typename Traits:: |
|
127 |
typedef typename Traits::Flow Flow; |
|
128 | 128 |
|
129 | 129 |
///The type of the flow map. |
130 | 130 |
typedef typename Traits::FlowMap FlowMap; |
131 | 131 |
///The type of the elevator. |
132 | 132 |
typedef typename Traits::Elevator Elevator; |
133 | 133 |
///The type of the tolerance. |
134 | 134 |
typedef typename Traits::Tolerance Tolerance; |
135 | 135 |
|
136 | 136 |
private: |
137 | 137 |
|
138 | 138 |
TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); |
139 | 139 |
|
140 | 140 |
const Digraph& _graph; |
141 | 141 |
const CapacityMap* _capacity; |
142 | 142 |
|
143 | 143 |
int _node_num; |
144 | 144 |
|
145 | 145 |
Node _source, _target; |
146 | 146 |
|
147 | 147 |
FlowMap* _flow; |
148 | 148 |
bool _local_flow; |
149 | 149 |
|
150 | 150 |
Elevator* _level; |
151 | 151 |
bool _local_level; |
152 | 152 |
|
153 |
typedef typename Digraph::template NodeMap< |
|
153 |
typedef typename Digraph::template NodeMap<Flow> ExcessMap; |
|
154 | 154 |
ExcessMap* _excess; |
155 | 155 |
|
156 | 156 |
Tolerance _tolerance; |
157 | 157 |
|
158 | 158 |
bool _phase; |
159 | 159 |
|
160 | 160 |
|
161 | 161 |
void createStructures() { |
162 | 162 |
_node_num = countNodes(_graph); |
163 | 163 |
|
164 | 164 |
if (!_flow) { |
165 | 165 |
_flow = Traits::createFlowMap(_graph); |
166 | 166 |
_local_flow = true; |
167 | 167 |
} |
168 | 168 |
if (!_level) { |
169 | 169 |
_level = Traits::createElevator(_graph, _node_num); |
170 | 170 |
_local_level = true; |
171 | 171 |
} |
172 | 172 |
if (!_excess) { |
173 | 173 |
_excess = new ExcessMap(_graph); |
174 | 174 |
} |
175 | 175 |
} |
176 | 176 |
|
177 | 177 |
void destroyStructures() { |
178 | 178 |
if (_local_flow) { |
179 | 179 |
delete _flow; |
180 | 180 |
} |
181 | 181 |
if (_local_level) { |
182 | 182 |
delete _level; |
183 | 183 |
} |
184 | 184 |
if (_excess) { |
185 | 185 |
delete _excess; |
186 | 186 |
} |
187 | 187 |
} |
188 | 188 |
|
189 | 189 |
public: |
190 | 190 |
|
191 | 191 |
typedef Preflow Create; |
192 | 192 |
|
193 | 193 |
///\name Named Template Parameters |
194 | 194 |
|
195 | 195 |
///@{ |
196 | 196 |
|
197 | 197 |
template <typename _FlowMap> |
198 | 198 |
struct SetFlowMapTraits : public Traits { |
199 | 199 |
typedef _FlowMap FlowMap; |
200 | 200 |
static FlowMap *createFlowMap(const Digraph&) { |
201 | 201 |
LEMON_ASSERT(false, "FlowMap is not initialized"); |
... | ... |
@@ -424,291 +424,291 @@ |
424 | 424 |
_level->initNewLevel(); |
425 | 425 |
std::vector<Node> nqueue; |
426 | 426 |
for (int i = 0; i < int(queue.size()); ++i) { |
427 | 427 |
Node n = queue[i]; |
428 | 428 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
429 | 429 |
Node u = _graph.source(e); |
430 | 430 |
if (!reached[u] && _tolerance.positive((*_capacity)[e])) { |
431 | 431 |
reached.set(u, true); |
432 | 432 |
_level->initAddItem(u); |
433 | 433 |
nqueue.push_back(u); |
434 | 434 |
} |
435 | 435 |
} |
436 | 436 |
} |
437 | 437 |
queue.swap(nqueue); |
438 | 438 |
} |
439 | 439 |
_level->initFinish(); |
440 | 440 |
|
441 | 441 |
for (OutArcIt e(_graph, _source); e != INVALID; ++e) { |
442 | 442 |
if (_tolerance.positive((*_capacity)[e])) { |
443 | 443 |
Node u = _graph.target(e); |
444 | 444 |
if ((*_level)[u] == _level->maxLevel()) continue; |
445 | 445 |
_flow->set(e, (*_capacity)[e]); |
446 | 446 |
_excess->set(u, (*_excess)[u] + (*_capacity)[e]); |
447 | 447 |
if (u != _target && !_level->active(u)) { |
448 | 448 |
_level->activate(u); |
449 | 449 |
} |
450 | 450 |
} |
451 | 451 |
} |
452 | 452 |
} |
453 | 453 |
|
454 | 454 |
/// \brief Initializes the internal data structures using the |
455 | 455 |
/// given flow map. |
456 | 456 |
/// |
457 | 457 |
/// Initializes the internal data structures and sets the initial |
458 | 458 |
/// flow to the given \c flowMap. The \c flowMap should contain a |
459 | 459 |
/// flow or at least a preflow, i.e. at each node excluding the |
460 | 460 |
/// source node the incoming flow should greater or equal to the |
461 | 461 |
/// outgoing flow. |
462 | 462 |
/// \return \c false if the given \c flowMap is not a preflow. |
463 | 463 |
template <typename FlowMap> |
464 | 464 |
bool init(const FlowMap& flowMap) { |
465 | 465 |
createStructures(); |
466 | 466 |
|
467 | 467 |
for (ArcIt e(_graph); e != INVALID; ++e) { |
468 | 468 |
_flow->set(e, flowMap[e]); |
469 | 469 |
} |
470 | 470 |
|
471 | 471 |
for (NodeIt n(_graph); n != INVALID; ++n) { |
472 |
|
|
472 |
Flow excess = 0; |
|
473 | 473 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
474 | 474 |
excess += (*_flow)[e]; |
475 | 475 |
} |
476 | 476 |
for (OutArcIt e(_graph, n); e != INVALID; ++e) { |
477 | 477 |
excess -= (*_flow)[e]; |
478 | 478 |
} |
479 | 479 |
if (excess < 0 && n != _source) return false; |
480 | 480 |
_excess->set(n, excess); |
481 | 481 |
} |
482 | 482 |
|
483 | 483 |
typename Digraph::template NodeMap<bool> reached(_graph, false); |
484 | 484 |
|
485 | 485 |
_level->initStart(); |
486 | 486 |
_level->initAddItem(_target); |
487 | 487 |
|
488 | 488 |
std::vector<Node> queue; |
489 | 489 |
reached.set(_source, true); |
490 | 490 |
|
491 | 491 |
queue.push_back(_target); |
492 | 492 |
reached.set(_target, true); |
493 | 493 |
while (!queue.empty()) { |
494 | 494 |
_level->initNewLevel(); |
495 | 495 |
std::vector<Node> nqueue; |
496 | 496 |
for (int i = 0; i < int(queue.size()); ++i) { |
497 | 497 |
Node n = queue[i]; |
498 | 498 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
499 | 499 |
Node u = _graph.source(e); |
500 | 500 |
if (!reached[u] && |
501 | 501 |
_tolerance.positive((*_capacity)[e] - (*_flow)[e])) { |
502 | 502 |
reached.set(u, true); |
503 | 503 |
_level->initAddItem(u); |
504 | 504 |
nqueue.push_back(u); |
505 | 505 |
} |
506 | 506 |
} |
507 | 507 |
for (OutArcIt e(_graph, n); e != INVALID; ++e) { |
508 | 508 |
Node v = _graph.target(e); |
509 | 509 |
if (!reached[v] && _tolerance.positive((*_flow)[e])) { |
510 | 510 |
reached.set(v, true); |
511 | 511 |
_level->initAddItem(v); |
512 | 512 |
nqueue.push_back(v); |
513 | 513 |
} |
514 | 514 |
} |
515 | 515 |
} |
516 | 516 |
queue.swap(nqueue); |
517 | 517 |
} |
518 | 518 |
_level->initFinish(); |
519 | 519 |
|
520 | 520 |
for (OutArcIt e(_graph, _source); e != INVALID; ++e) { |
521 |
|
|
521 |
Flow rem = (*_capacity)[e] - (*_flow)[e]; |
|
522 | 522 |
if (_tolerance.positive(rem)) { |
523 | 523 |
Node u = _graph.target(e); |
524 | 524 |
if ((*_level)[u] == _level->maxLevel()) continue; |
525 | 525 |
_flow->set(e, (*_capacity)[e]); |
526 | 526 |
_excess->set(u, (*_excess)[u] + rem); |
527 | 527 |
if (u != _target && !_level->active(u)) { |
528 | 528 |
_level->activate(u); |
529 | 529 |
} |
530 | 530 |
} |
531 | 531 |
} |
532 | 532 |
for (InArcIt e(_graph, _source); e != INVALID; ++e) { |
533 |
|
|
533 |
Flow rem = (*_flow)[e]; |
|
534 | 534 |
if (_tolerance.positive(rem)) { |
535 | 535 |
Node v = _graph.source(e); |
536 | 536 |
if ((*_level)[v] == _level->maxLevel()) continue; |
537 | 537 |
_flow->set(e, 0); |
538 | 538 |
_excess->set(v, (*_excess)[v] + rem); |
539 | 539 |
if (v != _target && !_level->active(v)) { |
540 | 540 |
_level->activate(v); |
541 | 541 |
} |
542 | 542 |
} |
543 | 543 |
} |
544 | 544 |
return true; |
545 | 545 |
} |
546 | 546 |
|
547 | 547 |
/// \brief Starts the first phase of the preflow algorithm. |
548 | 548 |
/// |
549 | 549 |
/// The preflow algorithm consists of two phases, this method runs |
550 | 550 |
/// the first phase. After the first phase the maximum flow value |
551 | 551 |
/// and a minimum value cut can already be computed, although a |
552 | 552 |
/// maximum flow is not yet obtained. So after calling this method |
553 | 553 |
/// \ref flowValue() returns the value of a maximum flow and \ref |
554 | 554 |
/// minCut() returns a minimum cut. |
555 | 555 |
/// \pre One of the \ref init() functions must be called before |
556 | 556 |
/// using this function. |
557 | 557 |
void startFirstPhase() { |
558 | 558 |
_phase = true; |
559 | 559 |
|
560 | 560 |
Node n = _level->highestActive(); |
561 | 561 |
int level = _level->highestActiveLevel(); |
562 | 562 |
while (n != INVALID) { |
563 | 563 |
int num = _node_num; |
564 | 564 |
|
565 | 565 |
while (num > 0 && n != INVALID) { |
566 |
|
|
566 |
Flow excess = (*_excess)[n]; |
|
567 | 567 |
int new_level = _level->maxLevel(); |
568 | 568 |
|
569 | 569 |
for (OutArcIt e(_graph, n); e != INVALID; ++e) { |
570 |
|
|
570 |
Flow rem = (*_capacity)[e] - (*_flow)[e]; |
|
571 | 571 |
if (!_tolerance.positive(rem)) continue; |
572 | 572 |
Node v = _graph.target(e); |
573 | 573 |
if ((*_level)[v] < level) { |
574 | 574 |
if (!_level->active(v) && v != _target) { |
575 | 575 |
_level->activate(v); |
576 | 576 |
} |
577 | 577 |
if (!_tolerance.less(rem, excess)) { |
578 | 578 |
_flow->set(e, (*_flow)[e] + excess); |
579 | 579 |
_excess->set(v, (*_excess)[v] + excess); |
580 | 580 |
excess = 0; |
581 | 581 |
goto no_more_push_1; |
582 | 582 |
} else { |
583 | 583 |
excess -= rem; |
584 | 584 |
_excess->set(v, (*_excess)[v] + rem); |
585 | 585 |
_flow->set(e, (*_capacity)[e]); |
586 | 586 |
} |
587 | 587 |
} else if (new_level > (*_level)[v]) { |
588 | 588 |
new_level = (*_level)[v]; |
589 | 589 |
} |
590 | 590 |
} |
591 | 591 |
|
592 | 592 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
593 |
|
|
593 |
Flow rem = (*_flow)[e]; |
|
594 | 594 |
if (!_tolerance.positive(rem)) continue; |
595 | 595 |
Node v = _graph.source(e); |
596 | 596 |
if ((*_level)[v] < level) { |
597 | 597 |
if (!_level->active(v) && v != _target) { |
598 | 598 |
_level->activate(v); |
599 | 599 |
} |
600 | 600 |
if (!_tolerance.less(rem, excess)) { |
601 | 601 |
_flow->set(e, (*_flow)[e] - excess); |
602 | 602 |
_excess->set(v, (*_excess)[v] + excess); |
603 | 603 |
excess = 0; |
604 | 604 |
goto no_more_push_1; |
605 | 605 |
} else { |
606 | 606 |
excess -= rem; |
607 | 607 |
_excess->set(v, (*_excess)[v] + rem); |
608 | 608 |
_flow->set(e, 0); |
609 | 609 |
} |
610 | 610 |
} else if (new_level > (*_level)[v]) { |
611 | 611 |
new_level = (*_level)[v]; |
612 | 612 |
} |
613 | 613 |
} |
614 | 614 |
|
615 | 615 |
no_more_push_1: |
616 | 616 |
|
617 | 617 |
_excess->set(n, excess); |
618 | 618 |
|
619 | 619 |
if (excess != 0) { |
620 | 620 |
if (new_level + 1 < _level->maxLevel()) { |
621 | 621 |
_level->liftHighestActive(new_level + 1); |
622 | 622 |
} else { |
623 | 623 |
_level->liftHighestActiveToTop(); |
624 | 624 |
} |
625 | 625 |
if (_level->emptyLevel(level)) { |
626 | 626 |
_level->liftToTop(level); |
627 | 627 |
} |
628 | 628 |
} else { |
629 | 629 |
_level->deactivate(n); |
630 | 630 |
} |
631 | 631 |
|
632 | 632 |
n = _level->highestActive(); |
633 | 633 |
level = _level->highestActiveLevel(); |
634 | 634 |
--num; |
635 | 635 |
} |
636 | 636 |
|
637 | 637 |
num = _node_num * 20; |
638 | 638 |
while (num > 0 && n != INVALID) { |
639 |
|
|
639 |
Flow excess = (*_excess)[n]; |
|
640 | 640 |
int new_level = _level->maxLevel(); |
641 | 641 |
|
642 | 642 |
for (OutArcIt e(_graph, n); e != INVALID; ++e) { |
643 |
|
|
643 |
Flow rem = (*_capacity)[e] - (*_flow)[e]; |
|
644 | 644 |
if (!_tolerance.positive(rem)) continue; |
645 | 645 |
Node v = _graph.target(e); |
646 | 646 |
if ((*_level)[v] < level) { |
647 | 647 |
if (!_level->active(v) && v != _target) { |
648 | 648 |
_level->activate(v); |
649 | 649 |
} |
650 | 650 |
if (!_tolerance.less(rem, excess)) { |
651 | 651 |
_flow->set(e, (*_flow)[e] + excess); |
652 | 652 |
_excess->set(v, (*_excess)[v] + excess); |
653 | 653 |
excess = 0; |
654 | 654 |
goto no_more_push_2; |
655 | 655 |
} else { |
656 | 656 |
excess -= rem; |
657 | 657 |
_excess->set(v, (*_excess)[v] + rem); |
658 | 658 |
_flow->set(e, (*_capacity)[e]); |
659 | 659 |
} |
660 | 660 |
} else if (new_level > (*_level)[v]) { |
661 | 661 |
new_level = (*_level)[v]; |
662 | 662 |
} |
663 | 663 |
} |
664 | 664 |
|
665 | 665 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
666 |
|
|
666 |
Flow rem = (*_flow)[e]; |
|
667 | 667 |
if (!_tolerance.positive(rem)) continue; |
668 | 668 |
Node v = _graph.source(e); |
669 | 669 |
if ((*_level)[v] < level) { |
670 | 670 |
if (!_level->active(v) && v != _target) { |
671 | 671 |
_level->activate(v); |
672 | 672 |
} |
673 | 673 |
if (!_tolerance.less(rem, excess)) { |
674 | 674 |
_flow->set(e, (*_flow)[e] - excess); |
675 | 675 |
_excess->set(v, (*_excess)[v] + excess); |
676 | 676 |
excess = 0; |
677 | 677 |
goto no_more_push_2; |
678 | 678 |
} else { |
679 | 679 |
excess -= rem; |
680 | 680 |
_excess->set(v, (*_excess)[v] + rem); |
681 | 681 |
_flow->set(e, 0); |
682 | 682 |
} |
683 | 683 |
} else if (new_level > (*_level)[v]) { |
684 | 684 |
new_level = (*_level)[v]; |
685 | 685 |
} |
686 | 686 |
} |
687 | 687 |
|
688 | 688 |
no_more_push_2: |
689 | 689 |
|
690 | 690 |
_excess->set(n, excess); |
691 | 691 |
|
692 | 692 |
if (excess != 0) { |
693 | 693 |
if (new_level + 1 < _level->maxLevel()) { |
694 | 694 |
_level->liftActiveOn(level, new_level + 1); |
695 | 695 |
} else { |
696 | 696 |
_level->liftActiveToTop(level); |
697 | 697 |
} |
698 | 698 |
if (_level->emptyLevel(level)) { |
699 | 699 |
_level->liftToTop(level); |
700 | 700 |
} |
701 | 701 |
} else { |
702 | 702 |
_level->deactivate(n); |
703 | 703 |
} |
704 | 704 |
|
705 | 705 |
while (level >= 0 && _level->activeFree(level)) { |
706 | 706 |
--level; |
707 | 707 |
} |
708 | 708 |
if (level == -1) { |
709 | 709 |
n = _level->highestActive(); |
710 | 710 |
level = _level->highestActiveLevel(); |
711 | 711 |
} else { |
712 | 712 |
n = _level->activeOn(level); |
713 | 713 |
} |
714 | 714 |
--num; |
... | ... |
@@ -732,227 +732,227 @@ |
732 | 732 |
for (NodeIt n(_graph); n != INVALID; ++n) { |
733 | 733 |
reached.set(n, (*_level)[n] < _level->maxLevel()); |
734 | 734 |
} |
735 | 735 |
|
736 | 736 |
_level->initStart(); |
737 | 737 |
_level->initAddItem(_source); |
738 | 738 |
|
739 | 739 |
std::vector<Node> queue; |
740 | 740 |
queue.push_back(_source); |
741 | 741 |
reached.set(_source, true); |
742 | 742 |
|
743 | 743 |
while (!queue.empty()) { |
744 | 744 |
_level->initNewLevel(); |
745 | 745 |
std::vector<Node> nqueue; |
746 | 746 |
for (int i = 0; i < int(queue.size()); ++i) { |
747 | 747 |
Node n = queue[i]; |
748 | 748 |
for (OutArcIt e(_graph, n); e != INVALID; ++e) { |
749 | 749 |
Node v = _graph.target(e); |
750 | 750 |
if (!reached[v] && _tolerance.positive((*_flow)[e])) { |
751 | 751 |
reached.set(v, true); |
752 | 752 |
_level->initAddItem(v); |
753 | 753 |
nqueue.push_back(v); |
754 | 754 |
} |
755 | 755 |
} |
756 | 756 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
757 | 757 |
Node u = _graph.source(e); |
758 | 758 |
if (!reached[u] && |
759 | 759 |
_tolerance.positive((*_capacity)[e] - (*_flow)[e])) { |
760 | 760 |
reached.set(u, true); |
761 | 761 |
_level->initAddItem(u); |
762 | 762 |
nqueue.push_back(u); |
763 | 763 |
} |
764 | 764 |
} |
765 | 765 |
} |
766 | 766 |
queue.swap(nqueue); |
767 | 767 |
} |
768 | 768 |
_level->initFinish(); |
769 | 769 |
|
770 | 770 |
for (NodeIt n(_graph); n != INVALID; ++n) { |
771 | 771 |
if (!reached[n]) { |
772 | 772 |
_level->dirtyTopButOne(n); |
773 | 773 |
} else if ((*_excess)[n] > 0 && _target != n) { |
774 | 774 |
_level->activate(n); |
775 | 775 |
} |
776 | 776 |
} |
777 | 777 |
|
778 | 778 |
Node n; |
779 | 779 |
while ((n = _level->highestActive()) != INVALID) { |
780 |
|
|
780 |
Flow excess = (*_excess)[n]; |
|
781 | 781 |
int level = _level->highestActiveLevel(); |
782 | 782 |
int new_level = _level->maxLevel(); |
783 | 783 |
|
784 | 784 |
for (OutArcIt e(_graph, n); e != INVALID; ++e) { |
785 |
|
|
785 |
Flow rem = (*_capacity)[e] - (*_flow)[e]; |
|
786 | 786 |
if (!_tolerance.positive(rem)) continue; |
787 | 787 |
Node v = _graph.target(e); |
788 | 788 |
if ((*_level)[v] < level) { |
789 | 789 |
if (!_level->active(v) && v != _source) { |
790 | 790 |
_level->activate(v); |
791 | 791 |
} |
792 | 792 |
if (!_tolerance.less(rem, excess)) { |
793 | 793 |
_flow->set(e, (*_flow)[e] + excess); |
794 | 794 |
_excess->set(v, (*_excess)[v] + excess); |
795 | 795 |
excess = 0; |
796 | 796 |
goto no_more_push; |
797 | 797 |
} else { |
798 | 798 |
excess -= rem; |
799 | 799 |
_excess->set(v, (*_excess)[v] + rem); |
800 | 800 |
_flow->set(e, (*_capacity)[e]); |
801 | 801 |
} |
802 | 802 |
} else if (new_level > (*_level)[v]) { |
803 | 803 |
new_level = (*_level)[v]; |
804 | 804 |
} |
805 | 805 |
} |
806 | 806 |
|
807 | 807 |
for (InArcIt e(_graph, n); e != INVALID; ++e) { |
808 |
|
|
808 |
Flow rem = (*_flow)[e]; |
|
809 | 809 |
if (!_tolerance.positive(rem)) continue; |
810 | 810 |
Node v = _graph.source(e); |
811 | 811 |
if ((*_level)[v] < level) { |
812 | 812 |
if (!_level->active(v) && v != _source) { |
813 | 813 |
_level->activate(v); |
814 | 814 |
} |
815 | 815 |
if (!_tolerance.less(rem, excess)) { |
816 | 816 |
_flow->set(e, (*_flow)[e] - excess); |
817 | 817 |
_excess->set(v, (*_excess)[v] + excess); |
818 | 818 |
excess = 0; |
819 | 819 |
goto no_more_push; |
820 | 820 |
} else { |
821 | 821 |
excess -= rem; |
822 | 822 |
_excess->set(v, (*_excess)[v] + rem); |
823 | 823 |
_flow->set(e, 0); |
824 | 824 |
} |
825 | 825 |
} else if (new_level > (*_level)[v]) { |
826 | 826 |
new_level = (*_level)[v]; |
827 | 827 |
} |
828 | 828 |
} |
829 | 829 |
|
830 | 830 |
no_more_push: |
831 | 831 |
|
832 | 832 |
_excess->set(n, excess); |
833 | 833 |
|
834 | 834 |
if (excess != 0) { |
835 | 835 |
if (new_level + 1 < _level->maxLevel()) { |
836 | 836 |
_level->liftHighestActive(new_level + 1); |
837 | 837 |
} else { |
838 | 838 |
// Calculation error |
839 | 839 |
_level->liftHighestActiveToTop(); |
840 | 840 |
} |
841 | 841 |
if (_level->emptyLevel(level)) { |
842 | 842 |
// Calculation error |
843 | 843 |
_level->liftToTop(level); |
844 | 844 |
} |
845 | 845 |
} else { |
846 | 846 |
_level->deactivate(n); |
847 | 847 |
} |
848 | 848 |
|
849 | 849 |
} |
850 | 850 |
} |
851 | 851 |
|
852 | 852 |
/// \brief Runs the preflow algorithm. |
853 | 853 |
/// |
854 | 854 |
/// Runs the preflow algorithm. |
855 | 855 |
/// \note pf.run() is just a shortcut of the following code. |
856 | 856 |
/// \code |
857 | 857 |
/// pf.init(); |
858 | 858 |
/// pf.startFirstPhase(); |
859 | 859 |
/// pf.startSecondPhase(); |
860 | 860 |
/// \endcode |
861 | 861 |
void run() { |
862 | 862 |
init(); |
863 | 863 |
startFirstPhase(); |
864 | 864 |
startSecondPhase(); |
865 | 865 |
} |
866 | 866 |
|
867 | 867 |
/// \brief Runs the preflow algorithm to compute the minimum cut. |
868 | 868 |
/// |
869 | 869 |
/// Runs the preflow algorithm to compute the minimum cut. |
870 | 870 |
/// \note pf.runMinCut() is just a shortcut of the following code. |
871 | 871 |
/// \code |
872 | 872 |
/// pf.init(); |
873 | 873 |
/// pf.startFirstPhase(); |
874 | 874 |
/// \endcode |
875 | 875 |
void runMinCut() { |
876 | 876 |
init(); |
877 | 877 |
startFirstPhase(); |
878 | 878 |
} |
879 | 879 |
|
880 | 880 |
/// @} |
881 | 881 |
|
882 | 882 |
/// \name Query Functions |
883 | 883 |
/// The results of the preflow algorithm can be obtained using these |
884 | 884 |
/// functions.\n |
885 | 885 |
/// Either one of the \ref run() "run*()" functions or one of the |
886 | 886 |
/// \ref startFirstPhase() "start*()" functions should be called |
887 | 887 |
/// before using them. |
888 | 888 |
|
889 | 889 |
///@{ |
890 | 890 |
|
891 | 891 |
/// \brief Returns the value of the maximum flow. |
892 | 892 |
/// |
893 | 893 |
/// Returns the value of the maximum flow by returning the excess |
894 | 894 |
/// of the target node. This value equals to the value of |
895 | 895 |
/// the maximum flow already after the first phase of the algorithm. |
896 | 896 |
/// |
897 | 897 |
/// \pre Either \ref run() or \ref init() must be called before |
898 | 898 |
/// using this function. |
899 |
|
|
899 |
Flow flowValue() const { |
|
900 | 900 |
return (*_excess)[_target]; |
901 | 901 |
} |
902 | 902 |
|
903 | 903 |
/// \brief Returns the flow on the given arc. |
904 | 904 |
/// |
905 | 905 |
/// Returns the flow on the given arc. This method can |
906 | 906 |
/// be called after the second phase of the algorithm. |
907 | 907 |
/// |
908 | 908 |
/// \pre Either \ref run() or \ref init() must be called before |
909 | 909 |
/// using this function. |
910 |
|
|
910 |
Flow flow(const Arc& arc) const { |
|
911 | 911 |
return (*_flow)[arc]; |
912 | 912 |
} |
913 | 913 |
|
914 | 914 |
/// \brief Returns a const reference to the flow map. |
915 | 915 |
/// |
916 | 916 |
/// Returns a const reference to the arc map storing the found flow. |
917 | 917 |
/// This method can be called after the second phase of the algorithm. |
918 | 918 |
/// |
919 | 919 |
/// \pre Either \ref run() or \ref init() must be called before |
920 | 920 |
/// using this function. |
921 | 921 |
const FlowMap& flowMap() const { |
922 | 922 |
return *_flow; |
923 | 923 |
} |
924 | 924 |
|
925 | 925 |
/// \brief Returns \c true when the node is on the source side of the |
926 | 926 |
/// minimum cut. |
927 | 927 |
/// |
928 | 928 |
/// Returns true when the node is on the source side of the found |
929 | 929 |
/// minimum cut. This method can be called both after running \ref |
930 | 930 |
/// startFirstPhase() and \ref startSecondPhase(). |
931 | 931 |
/// |
932 | 932 |
/// \pre Either \ref run() or \ref init() must be called before |
933 | 933 |
/// using this function. |
934 | 934 |
bool minCut(const Node& node) const { |
935 | 935 |
return ((*_level)[node] == _level->maxLevel()) == _phase; |
936 | 936 |
} |
937 | 937 |
|
938 | 938 |
/// \brief Gives back a minimum value cut. |
939 | 939 |
/// |
940 | 940 |
/// Sets \c cutMap to the characteristic vector of a minimum value |
941 | 941 |
/// cut. \c cutMap should be a \ref concepts::WriteMap "writable" |
942 | 942 |
/// node map with \c bool (or convertible) value type. |
943 | 943 |
/// |
944 | 944 |
/// This method can be called both after running \ref startFirstPhase() |
945 | 945 |
/// and \ref startSecondPhase(). The result after the second phase |
946 | 946 |
/// could be slightly different if inexact computation is used. |
947 | 947 |
/// |
948 | 948 |
/// \note This function calls \ref minCut() for each node, so it runs in |
949 | 949 |
/// \f$O(n)\f$ time. |
950 | 950 |
/// |
951 | 951 |
/// \pre Either \ref run() or \ref init() must be called before |
952 | 952 |
/// using this function. |
953 | 953 |
template <typename CutMap> |
954 | 954 |
void minCutMap(CutMap& cutMap) const { |
955 | 955 |
for (NodeIt n(_graph); n != INVALID; ++n) { |
956 | 956 |
cutMap.set(n, minCut(n)); |
957 | 957 |
} |
958 | 958 |
} |
... | ... |
@@ -12,120 +12,120 @@ |
12 | 12 |
* |
13 | 13 |
* This software is provided "AS IS" with no warranty of any kind, |
14 | 14 |
* express or implied, and with no claim as to its suitability for any |
15 | 15 |
* purpose. |
16 | 16 |
* |
17 | 17 |
*/ |
18 | 18 |
|
19 | 19 |
#include <iostream> |
20 | 20 |
|
21 | 21 |
#include "test_tools.h" |
22 | 22 |
#include <lemon/list_graph.h> |
23 | 23 |
#include <lemon/circulation.h> |
24 | 24 |
#include <lemon/lgf_reader.h> |
25 | 25 |
#include <lemon/concepts/digraph.h> |
26 | 26 |
#include <lemon/concepts/maps.h> |
27 | 27 |
|
28 | 28 |
using namespace lemon; |
29 | 29 |
|
30 | 30 |
char test_lgf[] = |
31 | 31 |
"@nodes\n" |
32 | 32 |
"label\n" |
33 | 33 |
"0\n" |
34 | 34 |
"1\n" |
35 | 35 |
"2\n" |
36 | 36 |
"3\n" |
37 | 37 |
"4\n" |
38 | 38 |
"5\n" |
39 | 39 |
"@arcs\n" |
40 | 40 |
" lcap ucap\n" |
41 | 41 |
"0 1 2 10\n" |
42 | 42 |
"0 2 2 6\n" |
43 | 43 |
"1 3 4 7\n" |
44 | 44 |
"1 4 0 5\n" |
45 | 45 |
"2 4 1 3\n" |
46 | 46 |
"3 5 3 8\n" |
47 | 47 |
"4 5 3 7\n" |
48 | 48 |
"@attributes\n" |
49 | 49 |
"source 0\n" |
50 | 50 |
"sink 5\n"; |
51 | 51 |
|
52 | 52 |
void checkCirculationCompile() |
53 | 53 |
{ |
54 | 54 |
typedef int VType; |
55 | 55 |
typedef concepts::Digraph Digraph; |
56 | 56 |
|
57 | 57 |
typedef Digraph::Node Node; |
58 | 58 |
typedef Digraph::Arc Arc; |
59 | 59 |
typedef concepts::ReadMap<Arc,VType> CapMap; |
60 |
typedef concepts::ReadMap<Node,VType> |
|
60 |
typedef concepts::ReadMap<Node,VType> SupplyMap; |
|
61 | 61 |
typedef concepts::ReadWriteMap<Arc,VType> FlowMap; |
62 | 62 |
typedef concepts::WriteMap<Node,bool> BarrierMap; |
63 | 63 |
|
64 | 64 |
typedef Elevator<Digraph, Digraph::Node> Elev; |
65 | 65 |
typedef LinkedElevator<Digraph, Digraph::Node> LinkedElev; |
66 | 66 |
|
67 | 67 |
Digraph g; |
68 | 68 |
Node n; |
69 | 69 |
Arc a; |
70 | 70 |
CapMap lcap, ucap; |
71 |
|
|
71 |
SupplyMap supply; |
|
72 | 72 |
FlowMap flow; |
73 | 73 |
BarrierMap bar; |
74 | 74 |
|
75 |
Circulation<Digraph, CapMap, CapMap, |
|
75 |
Circulation<Digraph, CapMap, CapMap, SupplyMap> |
|
76 | 76 |
::SetFlowMap<FlowMap> |
77 | 77 |
::SetElevator<Elev> |
78 | 78 |
::SetStandardElevator<LinkedElev> |
79 |
::Create circ_test(g,lcap,ucap, |
|
79 |
::Create circ_test(g,lcap,ucap,supply); |
|
80 | 80 |
|
81 |
circ_test.lowerCapMap(lcap); |
|
82 |
circ_test.upperCapMap(ucap); |
|
83 |
circ_test. |
|
81 |
circ_test.lowerMap(lcap); |
|
82 |
circ_test.upperMap(ucap); |
|
83 |
circ_test.supplyMap(supply); |
|
84 | 84 |
flow = circ_test.flowMap(); |
85 | 85 |
circ_test.flowMap(flow); |
86 | 86 |
|
87 | 87 |
circ_test.init(); |
88 | 88 |
circ_test.greedyInit(); |
89 | 89 |
circ_test.start(); |
90 | 90 |
circ_test.run(); |
91 | 91 |
|
92 | 92 |
circ_test.barrier(n); |
93 | 93 |
circ_test.barrierMap(bar); |
94 | 94 |
circ_test.flow(a); |
95 | 95 |
} |
96 | 96 |
|
97 | 97 |
template <class G, class LM, class UM, class DM> |
98 | 98 |
void checkCirculation(const G& g, const LM& lm, const UM& um, |
99 | 99 |
const DM& dm, bool find) |
100 | 100 |
{ |
101 | 101 |
Circulation<G, LM, UM, DM> circ(g, lm, um, dm); |
102 | 102 |
bool ret = circ.run(); |
103 | 103 |
if (find) { |
104 | 104 |
check(ret, "A feasible solution should have been found."); |
105 | 105 |
check(circ.checkFlow(), "The found flow is corrupt."); |
106 | 106 |
check(!circ.checkBarrier(), "A barrier should not have been found."); |
107 | 107 |
} else { |
108 | 108 |
check(!ret, "A feasible solution should not have been found."); |
109 | 109 |
check(circ.checkBarrier(), "The found barrier is corrupt."); |
110 | 110 |
} |
111 | 111 |
} |
112 | 112 |
|
113 | 113 |
int main (int, char*[]) |
114 | 114 |
{ |
115 | 115 |
typedef ListDigraph Digraph; |
116 | 116 |
DIGRAPH_TYPEDEFS(Digraph); |
117 | 117 |
|
118 | 118 |
Digraph g; |
119 | 119 |
IntArcMap lo(g), up(g); |
120 | 120 |
IntNodeMap delta(g, 0); |
121 | 121 |
Node s, t; |
122 | 122 |
|
123 | 123 |
std::istringstream input(test_lgf); |
124 | 124 |
DigraphReader<Digraph>(g,input). |
125 | 125 |
arcMap("lcap", lo). |
126 | 126 |
arcMap("ucap", up). |
127 | 127 |
node("source",s). |
128 | 128 |
node("sink",t). |
129 | 129 |
run(); |
130 | 130 |
|
131 | 131 |
delta[s] = 7; delta[t] = -7; |
0 comments (0 inline)