1.1 --- a/lemon/network_simplex.h Tue Mar 24 00:18:25 2009 +0100
1.2 +++ b/lemon/network_simplex.h Wed Mar 25 15:58:44 2009 +0100
1.3 @@ -22,7 +22,7 @@
1.4 /// \ingroup min_cost_flow
1.5 ///
1.6 /// \file
1.7 -/// \brief Network simplex algorithm for finding a minimum cost flow.
1.8 +/// \brief Network Simplex algorithm for finding a minimum cost flow.
1.9
1.10 #include <vector>
1.11 #include <limits>
1.12 @@ -36,80 +36,89 @@
1.13 /// \addtogroup min_cost_flow
1.14 /// @{
1.15
1.16 - /// \brief Implementation of the primal network simplex algorithm
1.17 + /// \brief Implementation of the primal Network Simplex algorithm
1.18 /// for finding a \ref min_cost_flow "minimum cost flow".
1.19 ///
1.20 - /// \ref NetworkSimplex implements the primal network simplex algorithm
1.21 + /// \ref NetworkSimplex implements the primal Network Simplex algorithm
1.22 /// for finding a \ref min_cost_flow "minimum cost flow".
1.23 ///
1.24 - /// \tparam Digraph The digraph type the algorithm runs on.
1.25 - /// \tparam LowerMap The type of the lower bound map.
1.26 - /// \tparam CapacityMap The type of the capacity (upper bound) map.
1.27 - /// \tparam CostMap The type of the cost (length) map.
1.28 - /// \tparam SupplyMap The type of the supply map.
1.29 + /// \tparam GR The digraph type the algorithm runs on.
1.30 + /// \tparam V The value type used in the algorithm.
1.31 + /// By default it is \c int.
1.32 ///
1.33 - /// \warning
1.34 - /// - Arc capacities and costs should be \e non-negative \e integers.
1.35 - /// - Supply values should be \e signed \e integers.
1.36 - /// - The value types of the maps should be convertible to each other.
1.37 - /// - \c CostMap::Value must be signed type.
1.38 + /// \warning \c V must be a signed integer type.
1.39 ///
1.40 - /// \note \ref NetworkSimplex provides five different pivot rule
1.41 - /// implementations that significantly affect the efficiency of the
1.42 - /// algorithm.
1.43 - /// By default "Block Search" pivot rule is used, which proved to be
1.44 - /// by far the most efficient according to our benchmark tests.
1.45 - /// However another pivot rule can be selected using \ref run()
1.46 - /// function with the proper parameter.
1.47 -#ifdef DOXYGEN
1.48 - template < typename Digraph,
1.49 - typename LowerMap,
1.50 - typename CapacityMap,
1.51 - typename CostMap,
1.52 - typename SupplyMap >
1.53 -
1.54 -#else
1.55 - template < typename Digraph,
1.56 - typename LowerMap = typename Digraph::template ArcMap<int>,
1.57 - typename CapacityMap = typename Digraph::template ArcMap<int>,
1.58 - typename CostMap = typename Digraph::template ArcMap<int>,
1.59 - typename SupplyMap = typename Digraph::template NodeMap<int> >
1.60 -#endif
1.61 + /// \note %NetworkSimplex provides five different pivot rule
1.62 + /// implementations. For more information see \ref PivotRule.
1.63 + template <typename GR, typename V = int>
1.64 class NetworkSimplex
1.65 {
1.66 - TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
1.67 + public:
1.68
1.69 - typedef typename CapacityMap::Value Capacity;
1.70 - typedef typename CostMap::Value Cost;
1.71 - typedef typename SupplyMap::Value Supply;
1.72 + /// The value type of the algorithm
1.73 + typedef V Value;
1.74 + /// The type of the flow map
1.75 + typedef typename GR::template ArcMap<Value> FlowMap;
1.76 + /// The type of the potential map
1.77 + typedef typename GR::template NodeMap<Value> PotentialMap;
1.78 +
1.79 + public:
1.80 +
1.81 + /// \brief Enum type for selecting the pivot rule.
1.82 + ///
1.83 + /// Enum type for selecting the pivot rule for the \ref run()
1.84 + /// function.
1.85 + ///
1.86 + /// \ref NetworkSimplex provides five different pivot rule
1.87 + /// implementations that significantly affect the running time
1.88 + /// of the algorithm.
1.89 + /// By default \ref BLOCK_SEARCH "Block Search" is used, which
1.90 + /// proved to be the most efficient and the most robust on various
1.91 + /// test inputs according to our benchmark tests.
1.92 + /// However another pivot rule can be selected using the \ref run()
1.93 + /// function with the proper parameter.
1.94 + enum PivotRule {
1.95 +
1.96 + /// The First Eligible pivot rule.
1.97 + /// The next eligible arc is selected in a wraparound fashion
1.98 + /// in every iteration.
1.99 + FIRST_ELIGIBLE,
1.100 +
1.101 + /// The Best Eligible pivot rule.
1.102 + /// The best eligible arc is selected in every iteration.
1.103 + BEST_ELIGIBLE,
1.104 +
1.105 + /// The Block Search pivot rule.
1.106 + /// A specified number of arcs are examined in every iteration
1.107 + /// in a wraparound fashion and the best eligible arc is selected
1.108 + /// from this block.
1.109 + BLOCK_SEARCH,
1.110 +
1.111 + /// The Candidate List pivot rule.
1.112 + /// In a major iteration a candidate list is built from eligible arcs
1.113 + /// in a wraparound fashion and in the following minor iterations
1.114 + /// the best eligible arc is selected from this list.
1.115 + CANDIDATE_LIST,
1.116 +
1.117 + /// The Altering Candidate List pivot rule.
1.118 + /// It is a modified version of the Candidate List method.
1.119 + /// It keeps only the several best eligible arcs from the former
1.120 + /// candidate list and extends this list in every iteration.
1.121 + ALTERING_LIST
1.122 + };
1.123 +
1.124 + private:
1.125 +
1.126 + TEMPLATE_DIGRAPH_TYPEDEFS(GR);
1.127 +
1.128 + typedef typename GR::template ArcMap<Value> ValueArcMap;
1.129 + typedef typename GR::template NodeMap<Value> ValueNodeMap;
1.130
1.131 typedef std::vector<Arc> ArcVector;
1.132 typedef std::vector<Node> NodeVector;
1.133 typedef std::vector<int> IntVector;
1.134 typedef std::vector<bool> BoolVector;
1.135 - typedef std::vector<Capacity> CapacityVector;
1.136 - typedef std::vector<Cost> CostVector;
1.137 - typedef std::vector<Supply> SupplyVector;
1.138 -
1.139 - public:
1.140 -
1.141 - /// The type of the flow map
1.142 - typedef typename Digraph::template ArcMap<Capacity> FlowMap;
1.143 - /// The type of the potential map
1.144 - typedef typename Digraph::template NodeMap<Cost> PotentialMap;
1.145 -
1.146 - public:
1.147 -
1.148 - /// Enum type for selecting the pivot rule used by \ref run()
1.149 - enum PivotRuleEnum {
1.150 - FIRST_ELIGIBLE_PIVOT,
1.151 - BEST_ELIGIBLE_PIVOT,
1.152 - BLOCK_SEARCH_PIVOT,
1.153 - CANDIDATE_LIST_PIVOT,
1.154 - ALTERING_LIST_PIVOT
1.155 - };
1.156 -
1.157 - private:
1.158 + typedef std::vector<Value> ValueVector;
1.159
1.160 // State constants for arcs
1.161 enum ArcStateEnum {
1.162 @@ -120,15 +129,19 @@
1.163
1.164 private:
1.165
1.166 - // References for the original data
1.167 - const Digraph &_graph;
1.168 - const LowerMap *_orig_lower;
1.169 - const CapacityMap &_orig_cap;
1.170 - const CostMap &_orig_cost;
1.171 - const SupplyMap *_orig_supply;
1.172 - Node _orig_source;
1.173 - Node _orig_target;
1.174 - Capacity _orig_flow_value;
1.175 + // Data related to the underlying digraph
1.176 + const GR &_graph;
1.177 + int _node_num;
1.178 + int _arc_num;
1.179 +
1.180 + // Parameters of the problem
1.181 + ValueArcMap *_plower;
1.182 + ValueArcMap *_pupper;
1.183 + ValueArcMap *_pcost;
1.184 + ValueNodeMap *_psupply;
1.185 + bool _pstsup;
1.186 + Node _psource, _ptarget;
1.187 + Value _pstflow;
1.188
1.189 // Result maps
1.190 FlowMap *_flow_map;
1.191 @@ -136,22 +149,18 @@
1.192 bool _local_flow;
1.193 bool _local_potential;
1.194
1.195 - // The number of nodes and arcs in the original graph
1.196 - int _node_num;
1.197 - int _arc_num;
1.198 -
1.199 - // Data structures for storing the graph
1.200 + // Data structures for storing the digraph
1.201 IntNodeMap _node_id;
1.202 ArcVector _arc_ref;
1.203 IntVector _source;
1.204 IntVector _target;
1.205
1.206 - // Node and arc maps
1.207 - CapacityVector _cap;
1.208 - CostVector _cost;
1.209 - CostVector _supply;
1.210 - CapacityVector _flow;
1.211 - CostVector _pi;
1.212 + // Node and arc data
1.213 + ValueVector _cap;
1.214 + ValueVector _cost;
1.215 + ValueVector _supply;
1.216 + ValueVector _flow;
1.217 + ValueVector _pi;
1.218
1.219 // Data for storing the spanning tree structure
1.220 IntVector _parent;
1.221 @@ -169,17 +178,11 @@
1.222 int in_arc, join, u_in, v_in, u_out, v_out;
1.223 int first, second, right, last;
1.224 int stem, par_stem, new_stem;
1.225 - Capacity delta;
1.226 + Value delta;
1.227
1.228 private:
1.229
1.230 - /// \brief Implementation of the "First Eligible" pivot rule for the
1.231 - /// \ref NetworkSimplex "network simplex" algorithm.
1.232 - ///
1.233 - /// This class implements the "First Eligible" pivot rule
1.234 - /// for the \ref NetworkSimplex "network simplex" algorithm.
1.235 - ///
1.236 - /// For more information see \ref NetworkSimplex::run().
1.237 + // Implementation of the First Eligible pivot rule
1.238 class FirstEligiblePivotRule
1.239 {
1.240 private:
1.241 @@ -187,9 +190,9 @@
1.242 // References to the NetworkSimplex class
1.243 const IntVector &_source;
1.244 const IntVector &_target;
1.245 - const CostVector &_cost;
1.246 + const ValueVector &_cost;
1.247 const IntVector &_state;
1.248 - const CostVector &_pi;
1.249 + const ValueVector &_pi;
1.250 int &_in_arc;
1.251 int _arc_num;
1.252
1.253 @@ -198,16 +201,16 @@
1.254
1.255 public:
1.256
1.257 - /// Constructor
1.258 + // Constructor
1.259 FirstEligiblePivotRule(NetworkSimplex &ns) :
1.260 _source(ns._source), _target(ns._target),
1.261 _cost(ns._cost), _state(ns._state), _pi(ns._pi),
1.262 _in_arc(ns.in_arc), _arc_num(ns._arc_num), _next_arc(0)
1.263 {}
1.264
1.265 - /// Find next entering arc
1.266 + // Find next entering arc
1.267 bool findEnteringArc() {
1.268 - Cost c;
1.269 + Value c;
1.270 for (int e = _next_arc; e < _arc_num; ++e) {
1.271 c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
1.272 if (c < 0) {
1.273 @@ -230,13 +233,7 @@
1.274 }; //class FirstEligiblePivotRule
1.275
1.276
1.277 - /// \brief Implementation of the "Best Eligible" pivot rule for the
1.278 - /// \ref NetworkSimplex "network simplex" algorithm.
1.279 - ///
1.280 - /// This class implements the "Best Eligible" pivot rule
1.281 - /// for the \ref NetworkSimplex "network simplex" algorithm.
1.282 - ///
1.283 - /// For more information see \ref NetworkSimplex::run().
1.284 + // Implementation of the Best Eligible pivot rule
1.285 class BestEligiblePivotRule
1.286 {
1.287 private:
1.288 @@ -244,24 +241,24 @@
1.289 // References to the NetworkSimplex class
1.290 const IntVector &_source;
1.291 const IntVector &_target;
1.292 - const CostVector &_cost;
1.293 + const ValueVector &_cost;
1.294 const IntVector &_state;
1.295 - const CostVector &_pi;
1.296 + const ValueVector &_pi;
1.297 int &_in_arc;
1.298 int _arc_num;
1.299
1.300 public:
1.301
1.302 - /// Constructor
1.303 + // Constructor
1.304 BestEligiblePivotRule(NetworkSimplex &ns) :
1.305 _source(ns._source), _target(ns._target),
1.306 _cost(ns._cost), _state(ns._state), _pi(ns._pi),
1.307 _in_arc(ns.in_arc), _arc_num(ns._arc_num)
1.308 {}
1.309
1.310 - /// Find next entering arc
1.311 + // Find next entering arc
1.312 bool findEnteringArc() {
1.313 - Cost c, min = 0;
1.314 + Value c, min = 0;
1.315 for (int e = 0; e < _arc_num; ++e) {
1.316 c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
1.317 if (c < min) {
1.318 @@ -275,13 +272,7 @@
1.319 }; //class BestEligiblePivotRule
1.320
1.321
1.322 - /// \brief Implementation of the "Block Search" pivot rule for the
1.323 - /// \ref NetworkSimplex "network simplex" algorithm.
1.324 - ///
1.325 - /// This class implements the "Block Search" pivot rule
1.326 - /// for the \ref NetworkSimplex "network simplex" algorithm.
1.327 - ///
1.328 - /// For more information see \ref NetworkSimplex::run().
1.329 + // Implementation of the Block Search pivot rule
1.330 class BlockSearchPivotRule
1.331 {
1.332 private:
1.333 @@ -289,9 +280,9 @@
1.334 // References to the NetworkSimplex class
1.335 const IntVector &_source;
1.336 const IntVector &_target;
1.337 - const CostVector &_cost;
1.338 + const ValueVector &_cost;
1.339 const IntVector &_state;
1.340 - const CostVector &_pi;
1.341 + const ValueVector &_pi;
1.342 int &_in_arc;
1.343 int _arc_num;
1.344
1.345 @@ -301,7 +292,7 @@
1.346
1.347 public:
1.348
1.349 - /// Constructor
1.350 + // Constructor
1.351 BlockSearchPivotRule(NetworkSimplex &ns) :
1.352 _source(ns._source), _target(ns._target),
1.353 _cost(ns._cost), _state(ns._state), _pi(ns._pi),
1.354 @@ -315,9 +306,9 @@
1.355 MIN_BLOCK_SIZE );
1.356 }
1.357
1.358 - /// Find next entering arc
1.359 + // Find next entering arc
1.360 bool findEnteringArc() {
1.361 - Cost c, min = 0;
1.362 + Value c, min = 0;
1.363 int cnt = _block_size;
1.364 int e, min_arc = _next_arc;
1.365 for (e = _next_arc; e < _arc_num; ++e) {
1.366 @@ -353,13 +344,7 @@
1.367 }; //class BlockSearchPivotRule
1.368
1.369
1.370 - /// \brief Implementation of the "Candidate List" pivot rule for the
1.371 - /// \ref NetworkSimplex "network simplex" algorithm.
1.372 - ///
1.373 - /// This class implements the "Candidate List" pivot rule
1.374 - /// for the \ref NetworkSimplex "network simplex" algorithm.
1.375 - ///
1.376 - /// For more information see \ref NetworkSimplex::run().
1.377 + // Implementation of the Candidate List pivot rule
1.378 class CandidateListPivotRule
1.379 {
1.380 private:
1.381 @@ -367,9 +352,9 @@
1.382 // References to the NetworkSimplex class
1.383 const IntVector &_source;
1.384 const IntVector &_target;
1.385 - const CostVector &_cost;
1.386 + const ValueVector &_cost;
1.387 const IntVector &_state;
1.388 - const CostVector &_pi;
1.389 + const ValueVector &_pi;
1.390 int &_in_arc;
1.391 int _arc_num;
1.392
1.393 @@ -403,7 +388,7 @@
1.394
1.395 /// Find next entering arc
1.396 bool findEnteringArc() {
1.397 - Cost min, c;
1.398 + Value min, c;
1.399 int e, min_arc = _next_arc;
1.400 if (_curr_length > 0 && _minor_count < _minor_limit) {
1.401 // Minor iteration: select the best eligible arc from the
1.402 @@ -464,13 +449,7 @@
1.403 }; //class CandidateListPivotRule
1.404
1.405
1.406 - /// \brief Implementation of the "Altering Candidate List" pivot rule
1.407 - /// for the \ref NetworkSimplex "network simplex" algorithm.
1.408 - ///
1.409 - /// This class implements the "Altering Candidate List" pivot rule
1.410 - /// for the \ref NetworkSimplex "network simplex" algorithm.
1.411 - ///
1.412 - /// For more information see \ref NetworkSimplex::run().
1.413 + // Implementation of the Altering Candidate List pivot rule
1.414 class AlteringListPivotRule
1.415 {
1.416 private:
1.417 @@ -478,9 +457,9 @@
1.418 // References to the NetworkSimplex class
1.419 const IntVector &_source;
1.420 const IntVector &_target;
1.421 - const CostVector &_cost;
1.422 + const ValueVector &_cost;
1.423 const IntVector &_state;
1.424 - const CostVector &_pi;
1.425 + const ValueVector &_pi;
1.426 int &_in_arc;
1.427 int _arc_num;
1.428
1.429 @@ -488,15 +467,15 @@
1.430 int _block_size, _head_length, _curr_length;
1.431 int _next_arc;
1.432 IntVector _candidates;
1.433 - CostVector _cand_cost;
1.434 + ValueVector _cand_cost;
1.435
1.436 // Functor class to compare arcs during sort of the candidate list
1.437 class SortFunc
1.438 {
1.439 private:
1.440 - const CostVector &_map;
1.441 + const ValueVector &_map;
1.442 public:
1.443 - SortFunc(const CostVector &map) : _map(map) {}
1.444 + SortFunc(const ValueVector &map) : _map(map) {}
1.445 bool operator()(int left, int right) {
1.446 return _map[left] > _map[right];
1.447 }
1.448 @@ -506,7 +485,7 @@
1.449
1.450 public:
1.451
1.452 - /// Constructor
1.453 + // Constructor
1.454 AlteringListPivotRule(NetworkSimplex &ns) :
1.455 _source(ns._source), _target(ns._target),
1.456 _cost(ns._cost), _state(ns._state), _pi(ns._pi),
1.457 @@ -527,7 +506,7 @@
1.458 _curr_length = 0;
1.459 }
1.460
1.461 - /// Find next entering arc
1.462 + // Find next entering arc
1.463 bool findEnteringArc() {
1.464 // Check the current candidate list
1.465 int e;
1.466 @@ -592,95 +571,23 @@
1.467
1.468 public:
1.469
1.470 - /// \brief General constructor (with lower bounds).
1.471 + /// \brief Constructor.
1.472 ///
1.473 - /// General constructor (with lower bounds).
1.474 + /// Constructor.
1.475 ///
1.476 /// \param graph The digraph the algorithm runs on.
1.477 - /// \param lower The lower bounds of the arcs.
1.478 - /// \param capacity The capacities (upper bounds) of the arcs.
1.479 - /// \param cost The cost (length) values of the arcs.
1.480 - /// \param supply The supply values of the nodes (signed).
1.481 - NetworkSimplex( const Digraph &graph,
1.482 - const LowerMap &lower,
1.483 - const CapacityMap &capacity,
1.484 - const CostMap &cost,
1.485 - const SupplyMap &supply ) :
1.486 - _graph(graph), _orig_lower(&lower), _orig_cap(capacity),
1.487 - _orig_cost(cost), _orig_supply(&supply),
1.488 + NetworkSimplex(const GR& graph) :
1.489 + _graph(graph),
1.490 + _plower(NULL), _pupper(NULL), _pcost(NULL),
1.491 + _psupply(NULL), _pstsup(false),
1.492 _flow_map(NULL), _potential_map(NULL),
1.493 _local_flow(false), _local_potential(false),
1.494 _node_id(graph)
1.495 - {}
1.496 -
1.497 - /// \brief General constructor (without lower bounds).
1.498 - ///
1.499 - /// General constructor (without lower bounds).
1.500 - ///
1.501 - /// \param graph The digraph the algorithm runs on.
1.502 - /// \param capacity The capacities (upper bounds) of the arcs.
1.503 - /// \param cost The cost (length) values of the arcs.
1.504 - /// \param supply The supply values of the nodes (signed).
1.505 - NetworkSimplex( const Digraph &graph,
1.506 - const CapacityMap &capacity,
1.507 - const CostMap &cost,
1.508 - const SupplyMap &supply ) :
1.509 - _graph(graph), _orig_lower(NULL), _orig_cap(capacity),
1.510 - _orig_cost(cost), _orig_supply(&supply),
1.511 - _flow_map(NULL), _potential_map(NULL),
1.512 - _local_flow(false), _local_potential(false),
1.513 - _node_id(graph)
1.514 - {}
1.515 -
1.516 - /// \brief Simple constructor (with lower bounds).
1.517 - ///
1.518 - /// Simple constructor (with lower bounds).
1.519 - ///
1.520 - /// \param graph The digraph the algorithm runs on.
1.521 - /// \param lower The lower bounds of the arcs.
1.522 - /// \param capacity The capacities (upper bounds) of the arcs.
1.523 - /// \param cost The cost (length) values of the arcs.
1.524 - /// \param s The source node.
1.525 - /// \param t The target node.
1.526 - /// \param flow_value The required amount of flow from node \c s
1.527 - /// to node \c t (i.e. the supply of \c s and the demand of \c t).
1.528 - NetworkSimplex( const Digraph &graph,
1.529 - const LowerMap &lower,
1.530 - const CapacityMap &capacity,
1.531 - const CostMap &cost,
1.532 - Node s, Node t,
1.533 - Capacity flow_value ) :
1.534 - _graph(graph), _orig_lower(&lower), _orig_cap(capacity),
1.535 - _orig_cost(cost), _orig_supply(NULL),
1.536 - _orig_source(s), _orig_target(t), _orig_flow_value(flow_value),
1.537 - _flow_map(NULL), _potential_map(NULL),
1.538 - _local_flow(false), _local_potential(false),
1.539 - _node_id(graph)
1.540 - {}
1.541 -
1.542 - /// \brief Simple constructor (without lower bounds).
1.543 - ///
1.544 - /// Simple constructor (without lower bounds).
1.545 - ///
1.546 - /// \param graph The digraph the algorithm runs on.
1.547 - /// \param capacity The capacities (upper bounds) of the arcs.
1.548 - /// \param cost The cost (length) values of the arcs.
1.549 - /// \param s The source node.
1.550 - /// \param t The target node.
1.551 - /// \param flow_value The required amount of flow from node \c s
1.552 - /// to node \c t (i.e. the supply of \c s and the demand of \c t).
1.553 - NetworkSimplex( const Digraph &graph,
1.554 - const CapacityMap &capacity,
1.555 - const CostMap &cost,
1.556 - Node s, Node t,
1.557 - Capacity flow_value ) :
1.558 - _graph(graph), _orig_lower(NULL), _orig_cap(capacity),
1.559 - _orig_cost(cost), _orig_supply(NULL),
1.560 - _orig_source(s), _orig_target(t), _orig_flow_value(flow_value),
1.561 - _flow_map(NULL), _potential_map(NULL),
1.562 - _local_flow(false), _local_potential(false),
1.563 - _node_id(graph)
1.564 - {}
1.565 + {
1.566 + LEMON_ASSERT(std::numeric_limits<Value>::is_integer &&
1.567 + std::numeric_limits<Value>::is_signed,
1.568 + "The value type of NetworkSimplex must be a signed integer");
1.569 + }
1.570
1.571 /// Destructor.
1.572 ~NetworkSimplex() {
1.573 @@ -688,12 +595,165 @@
1.574 if (_local_potential) delete _potential_map;
1.575 }
1.576
1.577 + /// \brief Set the lower bounds on the arcs.
1.578 + ///
1.579 + /// This function sets the lower bounds on the arcs.
1.580 + /// If neither this function nor \ref boundMaps() is used before
1.581 + /// calling \ref run(), the lower bounds will be set to zero
1.582 + /// on all arcs.
1.583 + ///
1.584 + /// \param map An arc map storing the lower bounds.
1.585 + /// Its \c Value type must be convertible to the \c Value type
1.586 + /// of the algorithm.
1.587 + ///
1.588 + /// \return <tt>(*this)</tt>
1.589 + template <typename LOWER>
1.590 + NetworkSimplex& lowerMap(const LOWER& map) {
1.591 + delete _plower;
1.592 + _plower = new ValueArcMap(_graph);
1.593 + for (ArcIt a(_graph); a != INVALID; ++a) {
1.594 + (*_plower)[a] = map[a];
1.595 + }
1.596 + return *this;
1.597 + }
1.598 +
1.599 + /// \brief Set the upper bounds (capacities) on the arcs.
1.600 + ///
1.601 + /// This function sets the upper bounds (capacities) on the arcs.
1.602 + /// If none of the functions \ref upperMap(), \ref capacityMap()
1.603 + /// and \ref boundMaps() is used before calling \ref run(),
1.604 + /// the upper bounds (capacities) will be set to
1.605 + /// \c std::numeric_limits<Value>::max() on all arcs.
1.606 + ///
1.607 + /// \param map An arc map storing the upper bounds.
1.608 + /// Its \c Value type must be convertible to the \c Value type
1.609 + /// of the algorithm.
1.610 + ///
1.611 + /// \return <tt>(*this)</tt>
1.612 + template<typename UPPER>
1.613 + NetworkSimplex& upperMap(const UPPER& map) {
1.614 + delete _pupper;
1.615 + _pupper = new ValueArcMap(_graph);
1.616 + for (ArcIt a(_graph); a != INVALID; ++a) {
1.617 + (*_pupper)[a] = map[a];
1.618 + }
1.619 + return *this;
1.620 + }
1.621 +
1.622 + /// \brief Set the upper bounds (capacities) on the arcs.
1.623 + ///
1.624 + /// This function sets the upper bounds (capacities) on the arcs.
1.625 + /// It is just an alias for \ref upperMap().
1.626 + ///
1.627 + /// \return <tt>(*this)</tt>
1.628 + template<typename CAP>
1.629 + NetworkSimplex& capacityMap(const CAP& map) {
1.630 + return upperMap(map);
1.631 + }
1.632 +
1.633 + /// \brief Set the lower and upper bounds on the arcs.
1.634 + ///
1.635 + /// This function sets the lower and upper bounds on the arcs.
1.636 + /// If neither this function nor \ref lowerMap() is used before
1.637 + /// calling \ref run(), the lower bounds will be set to zero
1.638 + /// on all arcs.
1.639 + /// If none of the functions \ref upperMap(), \ref capacityMap()
1.640 + /// and \ref boundMaps() is used before calling \ref run(),
1.641 + /// the upper bounds (capacities) will be set to
1.642 + /// \c std::numeric_limits<Value>::max() on all arcs.
1.643 + ///
1.644 + /// \param lower An arc map storing the lower bounds.
1.645 + /// \param upper An arc map storing the upper bounds.
1.646 + ///
1.647 + /// The \c Value type of the maps must be convertible to the
1.648 + /// \c Value type of the algorithm.
1.649 + ///
1.650 + /// \note This function is just a shortcut of calling \ref lowerMap()
1.651 + /// and \ref upperMap() separately.
1.652 + ///
1.653 + /// \return <tt>(*this)</tt>
1.654 + template <typename LOWER, typename UPPER>
1.655 + NetworkSimplex& boundMaps(const LOWER& lower, const UPPER& upper) {
1.656 + return lowerMap(lower).upperMap(upper);
1.657 + }
1.658 +
1.659 + /// \brief Set the costs of the arcs.
1.660 + ///
1.661 + /// This function sets the costs of the arcs.
1.662 + /// If it is not used before calling \ref run(), the costs
1.663 + /// will be set to \c 1 on all arcs.
1.664 + ///
1.665 + /// \param map An arc map storing the costs.
1.666 + /// Its \c Value type must be convertible to the \c Value type
1.667 + /// of the algorithm.
1.668 + ///
1.669 + /// \return <tt>(*this)</tt>
1.670 + template<typename COST>
1.671 + NetworkSimplex& costMap(const COST& map) {
1.672 + delete _pcost;
1.673 + _pcost = new ValueArcMap(_graph);
1.674 + for (ArcIt a(_graph); a != INVALID; ++a) {
1.675 + (*_pcost)[a] = map[a];
1.676 + }
1.677 + return *this;
1.678 + }
1.679 +
1.680 + /// \brief Set the supply values of the nodes.
1.681 + ///
1.682 + /// This function sets the supply values of the nodes.
1.683 + /// If neither this function nor \ref stSupply() is used before
1.684 + /// calling \ref run(), the supply of each node will be set to zero.
1.685 + /// (It makes sense only if non-zero lower bounds are given.)
1.686 + ///
1.687 + /// \param map A node map storing the supply values.
1.688 + /// Its \c Value type must be convertible to the \c Value type
1.689 + /// of the algorithm.
1.690 + ///
1.691 + /// \return <tt>(*this)</tt>
1.692 + template<typename SUP>
1.693 + NetworkSimplex& supplyMap(const SUP& map) {
1.694 + delete _psupply;
1.695 + _pstsup = false;
1.696 + _psupply = new ValueNodeMap(_graph);
1.697 + for (NodeIt n(_graph); n != INVALID; ++n) {
1.698 + (*_psupply)[n] = map[n];
1.699 + }
1.700 + return *this;
1.701 + }
1.702 +
1.703 + /// \brief Set single source and target nodes and a supply value.
1.704 + ///
1.705 + /// This function sets a single source node and a single target node
1.706 + /// and the required flow value.
1.707 + /// If neither this function nor \ref supplyMap() is used before
1.708 + /// calling \ref run(), the supply of each node will be set to zero.
1.709 + /// (It makes sense only if non-zero lower bounds are given.)
1.710 + ///
1.711 + /// \param s The source node.
1.712 + /// \param t The target node.
1.713 + /// \param k The required amount of flow from node \c s to node \c t
1.714 + /// (i.e. the supply of \c s and the demand of \c t).
1.715 + ///
1.716 + /// \return <tt>(*this)</tt>
1.717 + NetworkSimplex& stSupply(const Node& s, const Node& t, Value k) {
1.718 + delete _psupply;
1.719 + _psupply = NULL;
1.720 + _pstsup = true;
1.721 + _psource = s;
1.722 + _ptarget = t;
1.723 + _pstflow = k;
1.724 + return *this;
1.725 + }
1.726 +
1.727 /// \brief Set the flow map.
1.728 ///
1.729 /// This function sets the flow map.
1.730 + /// If it is not used before calling \ref run(), an instance will
1.731 + /// be allocated automatically. The destructor deallocates this
1.732 + /// automatically allocated map, of course.
1.733 ///
1.734 /// \return <tt>(*this)</tt>
1.735 - NetworkSimplex& flowMap(FlowMap &map) {
1.736 + NetworkSimplex& flowMap(FlowMap& map) {
1.737 if (_local_flow) {
1.738 delete _flow_map;
1.739 _local_flow = false;
1.740 @@ -704,10 +764,14 @@
1.741
1.742 /// \brief Set the potential map.
1.743 ///
1.744 - /// This function sets the potential map.
1.745 + /// This function sets the potential map, which is used for storing
1.746 + /// the dual solution.
1.747 + /// If it is not used before calling \ref run(), an instance will
1.748 + /// be allocated automatically. The destructor deallocates this
1.749 + /// automatically allocated map, of course.
1.750 ///
1.751 /// \return <tt>(*this)</tt>
1.752 - NetworkSimplex& potentialMap(PotentialMap &map) {
1.753 + NetworkSimplex& potentialMap(PotentialMap& map) {
1.754 if (_local_potential) {
1.755 delete _potential_map;
1.756 _local_potential = false;
1.757 @@ -716,47 +780,29 @@
1.758 return *this;
1.759 }
1.760
1.761 - /// \name Execution control
1.762 - /// The algorithm can be executed using the
1.763 - /// \ref lemon::NetworkSimplex::run() "run()" function.
1.764 + /// \name Execution Control
1.765 + /// The algorithm can be executed using \ref run().
1.766 +
1.767 /// @{
1.768
1.769 /// \brief Run the algorithm.
1.770 ///
1.771 /// This function runs the algorithm.
1.772 + /// The paramters can be specified using \ref lowerMap(),
1.773 + /// \ref upperMap(), \ref capacityMap(), \ref boundMaps(),
1.774 + /// \ref costMap(), \ref supplyMap() and \ref stSupply()
1.775 + /// functions. For example,
1.776 + /// \code
1.777 + /// NetworkSimplex<ListDigraph> ns(graph);
1.778 + /// ns.boundMaps(lower, upper).costMap(cost)
1.779 + /// .supplyMap(sup).run();
1.780 + /// \endcode
1.781 ///
1.782 - /// \param pivot_rule The pivot rule that is used during the
1.783 - /// algorithm.
1.784 - ///
1.785 - /// The available pivot rules:
1.786 - ///
1.787 - /// - FIRST_ELIGIBLE_PIVOT The next eligible arc is selected in
1.788 - /// a wraparound fashion in every iteration
1.789 - /// (\ref FirstEligiblePivotRule).
1.790 - ///
1.791 - /// - BEST_ELIGIBLE_PIVOT The best eligible arc is selected in
1.792 - /// every iteration (\ref BestEligiblePivotRule).
1.793 - ///
1.794 - /// - BLOCK_SEARCH_PIVOT A specified number of arcs are examined in
1.795 - /// every iteration in a wraparound fashion and the best eligible
1.796 - /// arc is selected from this block (\ref BlockSearchPivotRule).
1.797 - ///
1.798 - /// - CANDIDATE_LIST_PIVOT In a major iteration a candidate list is
1.799 - /// built from eligible arcs in a wraparound fashion and in the
1.800 - /// following minor iterations the best eligible arc is selected
1.801 - /// from this list (\ref CandidateListPivotRule).
1.802 - ///
1.803 - /// - ALTERING_LIST_PIVOT It is a modified version of the
1.804 - /// "Candidate List" pivot rule. It keeps only the several best
1.805 - /// eligible arcs from the former candidate list and extends this
1.806 - /// list in every iteration (\ref AlteringListPivotRule).
1.807 - ///
1.808 - /// According to our comprehensive benchmark tests the "Block Search"
1.809 - /// pivot rule proved to be the fastest and the most robust on
1.810 - /// various test inputs. Thus it is the default option.
1.811 + /// \param pivot_rule The pivot rule that will be used during the
1.812 + /// algorithm. For more information see \ref PivotRule.
1.813 ///
1.814 /// \return \c true if a feasible flow can be found.
1.815 - bool run(PivotRuleEnum pivot_rule = BLOCK_SEARCH_PIVOT) {
1.816 + bool run(PivotRule pivot_rule = BLOCK_SEARCH) {
1.817 return init() && start(pivot_rule);
1.818 }
1.819
1.820 @@ -765,10 +811,53 @@
1.821 /// \name Query Functions
1.822 /// The results of the algorithm can be obtained using these
1.823 /// functions.\n
1.824 - /// \ref lemon::NetworkSimplex::run() "run()" must be called before
1.825 - /// using them.
1.826 + /// The \ref run() function must be called before using them.
1.827 +
1.828 /// @{
1.829
1.830 + /// \brief Return the total cost of the found flow.
1.831 + ///
1.832 + /// This function returns the total cost of the found flow.
1.833 + /// The complexity of the function is \f$ O(e) \f$.
1.834 + ///
1.835 + /// \note The return type of the function can be specified as a
1.836 + /// template parameter. For example,
1.837 + /// \code
1.838 + /// ns.totalCost<double>();
1.839 + /// \endcode
1.840 + /// It is useful if the total cost cannot be stored in the \c Value
1.841 + /// type of the algorithm, which is the default return type of the
1.842 + /// function.
1.843 + ///
1.844 + /// \pre \ref run() must be called before using this function.
1.845 + template <typename Num>
1.846 + Num totalCost() const {
1.847 + Num c = 0;
1.848 + if (_pcost) {
1.849 + for (ArcIt e(_graph); e != INVALID; ++e)
1.850 + c += (*_flow_map)[e] * (*_pcost)[e];
1.851 + } else {
1.852 + for (ArcIt e(_graph); e != INVALID; ++e)
1.853 + c += (*_flow_map)[e];
1.854 + }
1.855 + return c;
1.856 + }
1.857 +
1.858 +#ifndef DOXYGEN
1.859 + Value totalCost() const {
1.860 + return totalCost<Value>();
1.861 + }
1.862 +#endif
1.863 +
1.864 + /// \brief Return the flow on the given arc.
1.865 + ///
1.866 + /// This function returns the flow on the given arc.
1.867 + ///
1.868 + /// \pre \ref run() must be called before using this function.
1.869 + Value flow(const Arc& a) const {
1.870 + return (*_flow_map)[a];
1.871 + }
1.872 +
1.873 /// \brief Return a const reference to the flow map.
1.874 ///
1.875 /// This function returns a const reference to an arc map storing
1.876 @@ -779,48 +868,28 @@
1.877 return *_flow_map;
1.878 }
1.879
1.880 + /// \brief Return the potential (dual value) of the given node.
1.881 + ///
1.882 + /// This function returns the potential (dual value) of the
1.883 + /// given node.
1.884 + ///
1.885 + /// \pre \ref run() must be called before using this function.
1.886 + Value potential(const Node& n) const {
1.887 + return (*_potential_map)[n];
1.888 + }
1.889 +
1.890 /// \brief Return a const reference to the potential map
1.891 /// (the dual solution).
1.892 ///
1.893 /// This function returns a const reference to a node map storing
1.894 - /// the found potentials (the dual solution).
1.895 + /// the found potentials, which form the dual solution of the
1.896 + /// \ref min_cost_flow "minimum cost flow" problem.
1.897 ///
1.898 /// \pre \ref run() must be called before using this function.
1.899 const PotentialMap& potentialMap() const {
1.900 return *_potential_map;
1.901 }
1.902
1.903 - /// \brief Return the flow on the given arc.
1.904 - ///
1.905 - /// This function returns the flow on the given arc.
1.906 - ///
1.907 - /// \pre \ref run() must be called before using this function.
1.908 - Capacity flow(const Arc& arc) const {
1.909 - return (*_flow_map)[arc];
1.910 - }
1.911 -
1.912 - /// \brief Return the potential of the given node.
1.913 - ///
1.914 - /// This function returns the potential of the given node.
1.915 - ///
1.916 - /// \pre \ref run() must be called before using this function.
1.917 - Cost potential(const Node& node) const {
1.918 - return (*_potential_map)[node];
1.919 - }
1.920 -
1.921 - /// \brief Return the total cost of the found flow.
1.922 - ///
1.923 - /// This function returns the total cost of the found flow.
1.924 - /// The complexity of the function is \f$ O(e) \f$.
1.925 - ///
1.926 - /// \pre \ref run() must be called before using this function.
1.927 - Cost totalCost() const {
1.928 - Cost c = 0;
1.929 - for (ArcIt e(_graph); e != INVALID; ++e)
1.930 - c += (*_flow_map)[e] * _orig_cost[e];
1.931 - return c;
1.932 - }
1.933 -
1.934 /// @}
1.935
1.936 private:
1.937 @@ -842,6 +911,7 @@
1.938 _arc_num = countArcs(_graph);
1.939 int all_node_num = _node_num + 1;
1.940 int all_arc_num = _arc_num + _node_num;
1.941 + if (_node_num == 0) return false;
1.942
1.943 _arc_ref.resize(_arc_num);
1.944 _source.resize(all_arc_num);
1.945 @@ -864,12 +934,17 @@
1.946
1.947 // Initialize node related data
1.948 bool valid_supply = true;
1.949 - if (_orig_supply) {
1.950 - Supply sum = 0;
1.951 + if (!_pstsup && !_psupply) {
1.952 + _pstsup = true;
1.953 + _psource = _ptarget = NodeIt(_graph);
1.954 + _pstflow = 0;
1.955 + }
1.956 + if (_psupply) {
1.957 + Value sum = 0;
1.958 int i = 0;
1.959 for (NodeIt n(_graph); n != INVALID; ++n, ++i) {
1.960 _node_id[n] = i;
1.961 - _supply[i] = (*_orig_supply)[n];
1.962 + _supply[i] = (*_psupply)[n];
1.963 sum += _supply[i];
1.964 }
1.965 valid_supply = (sum == 0);
1.966 @@ -879,8 +954,8 @@
1.967 _node_id[n] = i;
1.968 _supply[i] = 0;
1.969 }
1.970 - _supply[_node_id[_orig_source]] = _orig_flow_value;
1.971 - _supply[_node_id[_orig_target]] = -_orig_flow_value;
1.972 + _supply[_node_id[_psource]] = _pstflow;
1.973 + _supply[_node_id[_ptarget]] = -_pstflow;
1.974 }
1.975 if (!valid_supply) return false;
1.976
1.977 @@ -904,18 +979,41 @@
1.978 }
1.979
1.980 // Initialize arc maps
1.981 - for (int i = 0; i != _arc_num; ++i) {
1.982 - Arc e = _arc_ref[i];
1.983 - _source[i] = _node_id[_graph.source(e)];
1.984 - _target[i] = _node_id[_graph.target(e)];
1.985 - _cost[i] = _orig_cost[e];
1.986 - _cap[i] = _orig_cap[e];
1.987 + if (_pupper && _pcost) {
1.988 + for (int i = 0; i != _arc_num; ++i) {
1.989 + Arc e = _arc_ref[i];
1.990 + _source[i] = _node_id[_graph.source(e)];
1.991 + _target[i] = _node_id[_graph.target(e)];
1.992 + _cap[i] = (*_pupper)[e];
1.993 + _cost[i] = (*_pcost)[e];
1.994 + }
1.995 + } else {
1.996 + for (int i = 0; i != _arc_num; ++i) {
1.997 + Arc e = _arc_ref[i];
1.998 + _source[i] = _node_id[_graph.source(e)];
1.999 + _target[i] = _node_id[_graph.target(e)];
1.1000 + }
1.1001 + if (_pupper) {
1.1002 + for (int i = 0; i != _arc_num; ++i)
1.1003 + _cap[i] = (*_pupper)[_arc_ref[i]];
1.1004 + } else {
1.1005 + Value val = std::numeric_limits<Value>::max();
1.1006 + for (int i = 0; i != _arc_num; ++i)
1.1007 + _cap[i] = val;
1.1008 + }
1.1009 + if (_pcost) {
1.1010 + for (int i = 0; i != _arc_num; ++i)
1.1011 + _cost[i] = (*_pcost)[_arc_ref[i]];
1.1012 + } else {
1.1013 + for (int i = 0; i != _arc_num; ++i)
1.1014 + _cost[i] = 1;
1.1015 + }
1.1016 }
1.1017
1.1018 // Remove non-zero lower bounds
1.1019 - if (_orig_lower) {
1.1020 + if (_plower) {
1.1021 for (int i = 0; i != _arc_num; ++i) {
1.1022 - Capacity c = (*_orig_lower)[_arc_ref[i]];
1.1023 + Value c = (*_plower)[_arc_ref[i]];
1.1024 if (c != 0) {
1.1025 _cap[i] -= c;
1.1026 _supply[_source[i]] -= c;
1.1027 @@ -925,8 +1023,8 @@
1.1028 }
1.1029
1.1030 // Add artificial arcs and initialize the spanning tree data structure
1.1031 - Cost max_cost = std::numeric_limits<Cost>::max() / 4;
1.1032 - Capacity max_cap = std::numeric_limits<Capacity>::max();
1.1033 + Value max_cap = std::numeric_limits<Value>::max();
1.1034 + Value max_cost = std::numeric_limits<Value>::max() / 4;
1.1035 for (int u = 0, e = _arc_num; u != _node_num; ++u, ++e) {
1.1036 _thread[u] = u + 1;
1.1037 _rev_thread[u + 1] = u;
1.1038 @@ -979,7 +1077,7 @@
1.1039 }
1.1040 delta = _cap[in_arc];
1.1041 int result = 0;
1.1042 - Capacity d;
1.1043 + Value d;
1.1044 int e;
1.1045
1.1046 // Search the cycle along the path form the first node to the root
1.1047 @@ -1017,7 +1115,7 @@
1.1048 void changeFlow(bool change) {
1.1049 // Augment along the cycle
1.1050 if (delta > 0) {
1.1051 - Capacity val = _state[in_arc] * delta;
1.1052 + Value val = _state[in_arc] * delta;
1.1053 _flow[in_arc] += val;
1.1054 for (int u = _source[in_arc]; u != join; u = _parent[u]) {
1.1055 _flow[_pred[u]] += _forward[u] ? -val : val;
1.1056 @@ -1158,7 +1256,7 @@
1.1057
1.1058 // Update potentials
1.1059 void updatePotential() {
1.1060 - Cost sigma = _forward[u_in] ?
1.1061 + Value sigma = _forward[u_in] ?
1.1062 _pi[v_in] - _pi[u_in] - _cost[_pred[u_in]] :
1.1063 _pi[v_in] - _pi[u_in] + _cost[_pred[u_in]];
1.1064 if (_succ_num[u_in] > _node_num / 2) {
1.1065 @@ -1181,28 +1279,28 @@
1.1066 }
1.1067
1.1068 // Execute the algorithm
1.1069 - bool start(PivotRuleEnum pivot_rule) {
1.1070 + bool start(PivotRule pivot_rule) {
1.1071 // Select the pivot rule implementation
1.1072 switch (pivot_rule) {
1.1073 - case FIRST_ELIGIBLE_PIVOT:
1.1074 + case FIRST_ELIGIBLE:
1.1075 return start<FirstEligiblePivotRule>();
1.1076 - case BEST_ELIGIBLE_PIVOT:
1.1077 + case BEST_ELIGIBLE:
1.1078 return start<BestEligiblePivotRule>();
1.1079 - case BLOCK_SEARCH_PIVOT:
1.1080 + case BLOCK_SEARCH:
1.1081 return start<BlockSearchPivotRule>();
1.1082 - case CANDIDATE_LIST_PIVOT:
1.1083 + case CANDIDATE_LIST:
1.1084 return start<CandidateListPivotRule>();
1.1085 - case ALTERING_LIST_PIVOT:
1.1086 + case ALTERING_LIST:
1.1087 return start<AlteringListPivotRule>();
1.1088 }
1.1089 return false;
1.1090 }
1.1091
1.1092 - template<class PivotRuleImplementation>
1.1093 + template <typename PivotRuleImpl>
1.1094 bool start() {
1.1095 - PivotRuleImplementation pivot(*this);
1.1096 + PivotRuleImpl pivot(*this);
1.1097
1.1098 - // Execute the network simplex algorithm
1.1099 + // Execute the Network Simplex algorithm
1.1100 while (pivot.findEnteringArc()) {
1.1101 findJoinNode();
1.1102 bool change = findLeavingArc();
1.1103 @@ -1219,10 +1317,10 @@
1.1104 }
1.1105
1.1106 // Copy flow values to _flow_map
1.1107 - if (_orig_lower) {
1.1108 + if (_plower) {
1.1109 for (int i = 0; i != _arc_num; ++i) {
1.1110 Arc e = _arc_ref[i];
1.1111 - _flow_map->set(e, (*_orig_lower)[e] + _flow[i]);
1.1112 + _flow_map->set(e, (*_plower)[e] + _flow[i]);
1.1113 }
1.1114 } else {
1.1115 for (int i = 0; i != _arc_num; ++i) {
2.1 --- a/test/min_cost_flow_test.cc Tue Mar 24 00:18:25 2009 +0100
2.2 +++ b/test/min_cost_flow_test.cc Wed Mar 25 15:58:44 2009 +0100
2.3 @@ -20,15 +20,9 @@
2.4 #include <fstream>
2.5
2.6 #include <lemon/list_graph.h>
2.7 -#include <lemon/smart_graph.h>
2.8 #include <lemon/lgf_reader.h>
2.9
2.10 -//#include <lemon/cycle_canceling.h>
2.11 -//#include <lemon/capacity_scaling.h>
2.12 -//#include <lemon/cost_scaling.h>
2.13 #include <lemon/network_simplex.h>
2.14 -//#include <lemon/min_cost_flow.h>
2.15 -//#include <lemon/min_cost_max_flow.h>
2.16
2.17 #include <lemon/concepts/digraph.h>
2.18 #include <lemon/concept_check.h>
2.19 @@ -93,36 +87,30 @@
2.20 void constraints() {
2.21 checkConcept<concepts::Digraph, GR>();
2.22
2.23 - MCF mcf_test1(g, lower, upper, cost, sup);
2.24 - MCF mcf_test2(g, upper, cost, sup);
2.25 - MCF mcf_test3(g, lower, upper, cost, n, n, k);
2.26 - MCF mcf_test4(g, upper, cost, n, n, k);
2.27 + MCF mcf(g);
2.28
2.29 - // TODO: This part should be enabled and the next part
2.30 - // should be removed if map copying is supported
2.31 -/*
2.32 - flow = mcf_test1.flowMap();
2.33 - mcf_test1.flowMap(flow);
2.34 + b = mcf.lowerMap(lower)
2.35 + .upperMap(upper)
2.36 + .capacityMap(upper)
2.37 + .boundMaps(lower, upper)
2.38 + .costMap(cost)
2.39 + .supplyMap(sup)
2.40 + .stSupply(n, n, k)
2.41 + .run();
2.42
2.43 - pot = mcf_test1.potentialMap();
2.44 - mcf_test1.potentialMap(pot);
2.45 -*/
2.46 -/**/
2.47 - const typename MCF::FlowMap &fm =
2.48 - mcf_test1.flowMap();
2.49 - mcf_test1.flowMap(flow);
2.50 - const typename MCF::PotentialMap &pm =
2.51 - mcf_test1.potentialMap();
2.52 - mcf_test1.potentialMap(pot);
2.53 + const typename MCF::FlowMap &fm = mcf.flowMap();
2.54 + const typename MCF::PotentialMap &pm = mcf.potentialMap();
2.55 +
2.56 + v = mcf.totalCost();
2.57 + double x = mcf.template totalCost<double>();
2.58 + v = mcf.flow(a);
2.59 + v = mcf.potential(n);
2.60 + mcf.flowMap(flow);
2.61 + mcf.potentialMap(pot);
2.62 +
2.63 ignore_unused_variable_warning(fm);
2.64 ignore_unused_variable_warning(pm);
2.65 -/**/
2.66 -
2.67 - mcf_test1.run();
2.68 -
2.69 - v = mcf_test1.totalCost();
2.70 - v = mcf_test1.flow(a);
2.71 - v = mcf_test1.potential(n);
2.72 + ignore_unused_variable_warning(x);
2.73 }
2.74
2.75 typedef typename GR::Node Node;
2.76 @@ -139,6 +127,7 @@
2.77 const Arc &a;
2.78 const Value &k;
2.79 Value v;
2.80 + bool b;
2.81
2.82 typename MCF::FlowMap &flow;
2.83 typename MCF::PotentialMap &pot;
2.84 @@ -172,7 +161,7 @@
2.85 }
2.86
2.87 // Check the feasibility of the given potentials (dual soluiton)
2.88 -// using the Complementary Slackness optimality condition
2.89 +// using the "Complementary Slackness" optimality condition
2.90 template < typename GR, typename LM, typename UM,
2.91 typename CM, typename FM, typename PM >
2.92 bool checkPotential( const GR& gr, const LM& lower, const UM& upper,
2.93 @@ -217,23 +206,14 @@
2.94 // Check the interfaces
2.95 {
2.96 typedef int Value;
2.97 - // This typedef should be enabled if the standard maps are
2.98 - // reference maps in the graph concepts
2.99 + // TODO: This typedef should be enabled if the standard maps are
2.100 + // reference maps in the graph concepts (See #190).
2.101 +/**/
2.102 //typedef concepts::Digraph GR;
2.103 typedef ListDigraph GR;
2.104 - typedef concepts::ReadMap<GR::Node, Value> NM;
2.105 - typedef concepts::ReadMap<GR::Arc, Value> AM;
2.106 -
2.107 - //checkConcept< McfClassConcept<GR, Value>,
2.108 - // CycleCanceling<GR, AM, AM, AM, NM> >();
2.109 - //checkConcept< McfClassConcept<GR, Value>,
2.110 - // CapacityScaling<GR, AM, AM, AM, NM> >();
2.111 - //checkConcept< McfClassConcept<GR, Value>,
2.112 - // CostScaling<GR, AM, AM, AM, NM> >();
2.113 +/**/
2.114 checkConcept< McfClassConcept<GR, Value>,
2.115 - NetworkSimplex<GR, AM, AM, AM, NM> >();
2.116 - //checkConcept< MinCostFlow<GR, Value>,
2.117 - // NetworkSimplex<GR, AM, AM, AM, NM> >();
2.118 + NetworkSimplex<GR, Value> >();
2.119 }
2.120
2.121 // Run various MCF tests
2.122 @@ -244,6 +224,7 @@
2.123 Digraph gr;
2.124 Digraph::ArcMap<int> c(gr), l1(gr), l2(gr), u(gr);
2.125 Digraph::NodeMap<int> s1(gr), s2(gr), s3(gr);
2.126 + ConstMap<Arc, int> cc(1), cu(std::numeric_limits<int>::max());
2.127 Node v, w;
2.128
2.129 std::istringstream input(test_lgf);
2.130 @@ -259,197 +240,50 @@
2.131 .node("target", w)
2.132 .run();
2.133
2.134 -/*
2.135 - // A. Test CapacityScaling with scaling
2.136 + // A. Test NetworkSimplex with the default pivot rule
2.137 {
2.138 - CapacityScaling<Digraph> mcf1(gr, u, c, s1);
2.139 - CapacityScaling<Digraph> mcf2(gr, u, c, v, w, 27);
2.140 - CapacityScaling<Digraph> mcf3(gr, u, c, s3);
2.141 - CapacityScaling<Digraph> mcf4(gr, l2, u, c, s1);
2.142 - CapacityScaling<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.143 - CapacityScaling<Digraph> mcf6(gr, l2, u, c, s3);
2.144 + NetworkSimplex<Digraph> mcf1(gr), mcf2(gr), mcf3(gr), mcf4(gr),
2.145 + mcf5(gr), mcf6(gr), mcf7(gr), mcf8(gr);
2.146
2.147 - checkMcf(mcf1, mcf1.run(), gr, l1, u, c, s1, true, 5240, "#A1");
2.148 - checkMcf(mcf2, mcf2.run(), gr, l1, u, c, s2, true, 7620, "#A2");
2.149 - checkMcf(mcf3, mcf3.run(), gr, l1, u, c, s3, true, 0, "#A3");
2.150 - checkMcf(mcf4, mcf4.run(), gr, l2, u, c, s1, true, 5970, "#A4");
2.151 - checkMcf(mcf5, mcf5.run(), gr, l2, u, c, s2, true, 8010, "#A5");
2.152 - checkMcf(mcf6, mcf6.run(), gr, l2, u, c, s3, false, 0, "#A6");
2.153 + checkMcf(mcf1, mcf1.upperMap(u).costMap(c).supplyMap(s1).run(),
2.154 + gr, l1, u, c, s1, true, 5240, "#A1");
2.155 + checkMcf(mcf2, mcf2.upperMap(u).costMap(c).stSupply(v, w, 27).run(),
2.156 + gr, l1, u, c, s2, true, 7620, "#A2");
2.157 + checkMcf(mcf3, mcf3.boundMaps(l2, u).costMap(c).supplyMap(s1).run(),
2.158 + gr, l2, u, c, s1, true, 5970, "#A3");
2.159 + checkMcf(mcf4, mcf4.boundMaps(l2, u).costMap(c).stSupply(v, w, 27).run(),
2.160 + gr, l2, u, c, s2, true, 8010, "#A4");
2.161 + checkMcf(mcf5, mcf5.supplyMap(s1).run(),
2.162 + gr, l1, cu, cc, s1, true, 74, "#A5");
2.163 + checkMcf(mcf6, mcf6.stSupply(v, w, 27).lowerMap(l2).run(),
2.164 + gr, l2, cu, cc, s2, true, 94, "#A6");
2.165 + checkMcf(mcf7, mcf7.run(),
2.166 + gr, l1, cu, cc, s3, true, 0, "#A7");
2.167 + checkMcf(mcf8, mcf8.boundMaps(l2, u).run(),
2.168 + gr, l2, u, cc, s3, false, 0, "#A8");
2.169 }
2.170
2.171 - // B. Test CapacityScaling without scaling
2.172 + // B. Test NetworkSimplex with each pivot rule
2.173 {
2.174 - CapacityScaling<Digraph> mcf1(gr, u, c, s1);
2.175 - CapacityScaling<Digraph> mcf2(gr, u, c, v, w, 27);
2.176 - CapacityScaling<Digraph> mcf3(gr, u, c, s3);
2.177 - CapacityScaling<Digraph> mcf4(gr, l2, u, c, s1);
2.178 - CapacityScaling<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.179 - CapacityScaling<Digraph> mcf6(gr, l2, u, c, s3);
2.180 + NetworkSimplex<Digraph> mcf1(gr), mcf2(gr), mcf3(gr), mcf4(gr), mcf5(gr);
2.181 + NetworkSimplex<Digraph>::PivotRule pr;
2.182
2.183 - checkMcf(mcf1, mcf1.run(false), gr, l1, u, c, s1, true, 5240, "#B1");
2.184 - checkMcf(mcf2, mcf2.run(false), gr, l1, u, c, s2, true, 7620, "#B2");
2.185 - checkMcf(mcf3, mcf3.run(false), gr, l1, u, c, s3, true, 0, "#B3");
2.186 - checkMcf(mcf4, mcf4.run(false), gr, l2, u, c, s1, true, 5970, "#B4");
2.187 - checkMcf(mcf5, mcf5.run(false), gr, l2, u, c, s2, true, 8010, "#B5");
2.188 - checkMcf(mcf6, mcf6.run(false), gr, l2, u, c, s3, false, 0, "#B6");
2.189 + pr = NetworkSimplex<Digraph>::FIRST_ELIGIBLE;
2.190 + checkMcf(mcf1, mcf1.boundMaps(l2, u).costMap(c).supplyMap(s1).run(pr),
2.191 + gr, l2, u, c, s1, true, 5970, "#B1");
2.192 + pr = NetworkSimplex<Digraph>::BEST_ELIGIBLE;
2.193 + checkMcf(mcf2, mcf2.boundMaps(l2, u).costMap(c).supplyMap(s1).run(pr),
2.194 + gr, l2, u, c, s1, true, 5970, "#B2");
2.195 + pr = NetworkSimplex<Digraph>::BLOCK_SEARCH;
2.196 + checkMcf(mcf3, mcf3.boundMaps(l2, u).costMap(c).supplyMap(s1).run(pr),
2.197 + gr, l2, u, c, s1, true, 5970, "#B3");
2.198 + pr = NetworkSimplex<Digraph>::CANDIDATE_LIST;
2.199 + checkMcf(mcf4, mcf4.boundMaps(l2, u).costMap(c).supplyMap(s1).run(pr),
2.200 + gr, l2, u, c, s1, true, 5970, "#B4");
2.201 + pr = NetworkSimplex<Digraph>::ALTERING_LIST;
2.202 + checkMcf(mcf5, mcf5.boundMaps(l2, u).costMap(c).supplyMap(s1).run(pr),
2.203 + gr, l2, u, c, s1, true, 5970, "#B5");
2.204 }
2.205
2.206 - // C. Test CostScaling using partial augment-relabel method
2.207 - {
2.208 - CostScaling<Digraph> mcf1(gr, u, c, s1);
2.209 - CostScaling<Digraph> mcf2(gr, u, c, v, w, 27);
2.210 - CostScaling<Digraph> mcf3(gr, u, c, s3);
2.211 - CostScaling<Digraph> mcf4(gr, l2, u, c, s1);
2.212 - CostScaling<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.213 - CostScaling<Digraph> mcf6(gr, l2, u, c, s3);
2.214 -
2.215 - checkMcf(mcf1, mcf1.run(), gr, l1, u, c, s1, true, 5240, "#C1");
2.216 - checkMcf(mcf2, mcf2.run(), gr, l1, u, c, s2, true, 7620, "#C2");
2.217 - checkMcf(mcf3, mcf3.run(), gr, l1, u, c, s3, true, 0, "#C3");
2.218 - checkMcf(mcf4, mcf4.run(), gr, l2, u, c, s1, true, 5970, "#C4");
2.219 - checkMcf(mcf5, mcf5.run(), gr, l2, u, c, s2, true, 8010, "#C5");
2.220 - checkMcf(mcf6, mcf6.run(), gr, l2, u, c, s3, false, 0, "#C6");
2.221 - }
2.222 -
2.223 - // D. Test CostScaling using push-relabel method
2.224 - {
2.225 - CostScaling<Digraph> mcf1(gr, u, c, s1);
2.226 - CostScaling<Digraph> mcf2(gr, u, c, v, w, 27);
2.227 - CostScaling<Digraph> mcf3(gr, u, c, s3);
2.228 - CostScaling<Digraph> mcf4(gr, l2, u, c, s1);
2.229 - CostScaling<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.230 - CostScaling<Digraph> mcf6(gr, l2, u, c, s3);
2.231 -
2.232 - checkMcf(mcf1, mcf1.run(false), gr, l1, u, c, s1, true, 5240, "#D1");
2.233 - checkMcf(mcf2, mcf2.run(false), gr, l1, u, c, s2, true, 7620, "#D2");
2.234 - checkMcf(mcf3, mcf3.run(false), gr, l1, u, c, s3, true, 0, "#D3");
2.235 - checkMcf(mcf4, mcf4.run(false), gr, l2, u, c, s1, true, 5970, "#D4");
2.236 - checkMcf(mcf5, mcf5.run(false), gr, l2, u, c, s2, true, 8010, "#D5");
2.237 - checkMcf(mcf6, mcf6.run(false), gr, l2, u, c, s3, false, 0, "#D6");
2.238 - }
2.239 -*/
2.240 -
2.241 - // E. Test NetworkSimplex with FIRST_ELIGIBLE_PIVOT
2.242 - {
2.243 - NetworkSimplex<Digraph>::PivotRuleEnum pr =
2.244 - NetworkSimplex<Digraph>::FIRST_ELIGIBLE_PIVOT;
2.245 - NetworkSimplex<Digraph> mcf1(gr, u, c, s1);
2.246 - NetworkSimplex<Digraph> mcf2(gr, u, c, v, w, 27);
2.247 - NetworkSimplex<Digraph> mcf3(gr, u, c, s3);
2.248 - NetworkSimplex<Digraph> mcf4(gr, l2, u, c, s1);
2.249 - NetworkSimplex<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.250 - NetworkSimplex<Digraph> mcf6(gr, l2, u, c, s3);
2.251 -
2.252 - checkMcf(mcf1, mcf1.run(pr), gr, l1, u, c, s1, true, 5240, "#E1");
2.253 - checkMcf(mcf2, mcf2.run(pr), gr, l1, u, c, s2, true, 7620, "#E2");
2.254 - checkMcf(mcf3, mcf3.run(pr), gr, l1, u, c, s3, true, 0, "#E3");
2.255 - checkMcf(mcf4, mcf4.run(pr), gr, l2, u, c, s1, true, 5970, "#E4");
2.256 - checkMcf(mcf5, mcf5.run(pr), gr, l2, u, c, s2, true, 8010, "#E5");
2.257 - checkMcf(mcf6, mcf6.run(pr), gr, l2, u, c, s3, false, 0, "#E6");
2.258 - }
2.259 -
2.260 - // F. Test NetworkSimplex with BEST_ELIGIBLE_PIVOT
2.261 - {
2.262 - NetworkSimplex<Digraph>::PivotRuleEnum pr =
2.263 - NetworkSimplex<Digraph>::BEST_ELIGIBLE_PIVOT;
2.264 - NetworkSimplex<Digraph> mcf1(gr, u, c, s1);
2.265 - NetworkSimplex<Digraph> mcf2(gr, u, c, v, w, 27);
2.266 - NetworkSimplex<Digraph> mcf3(gr, u, c, s3);
2.267 - NetworkSimplex<Digraph> mcf4(gr, l2, u, c, s1);
2.268 - NetworkSimplex<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.269 - NetworkSimplex<Digraph> mcf6(gr, l2, u, c, s3);
2.270 -
2.271 - checkMcf(mcf1, mcf1.run(pr), gr, l1, u, c, s1, true, 5240, "#F1");
2.272 - checkMcf(mcf2, mcf2.run(pr), gr, l1, u, c, s2, true, 7620, "#F2");
2.273 - checkMcf(mcf3, mcf3.run(pr), gr, l1, u, c, s3, true, 0, "#F3");
2.274 - checkMcf(mcf4, mcf4.run(pr), gr, l2, u, c, s1, true, 5970, "#F4");
2.275 - checkMcf(mcf5, mcf5.run(pr), gr, l2, u, c, s2, true, 8010, "#F5");
2.276 - checkMcf(mcf6, mcf6.run(pr), gr, l2, u, c, s3, false, 0, "#F6");
2.277 - }
2.278 -
2.279 - // G. Test NetworkSimplex with BLOCK_SEARCH_PIVOT
2.280 - {
2.281 - NetworkSimplex<Digraph>::PivotRuleEnum pr =
2.282 - NetworkSimplex<Digraph>::BLOCK_SEARCH_PIVOT;
2.283 - NetworkSimplex<Digraph> mcf1(gr, u, c, s1);
2.284 - NetworkSimplex<Digraph> mcf2(gr, u, c, v, w, 27);
2.285 - NetworkSimplex<Digraph> mcf3(gr, u, c, s3);
2.286 - NetworkSimplex<Digraph> mcf4(gr, l2, u, c, s1);
2.287 - NetworkSimplex<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.288 - NetworkSimplex<Digraph> mcf6(gr, l2, u, c, s3);
2.289 -
2.290 - checkMcf(mcf1, mcf1.run(pr), gr, l1, u, c, s1, true, 5240, "#G1");
2.291 - checkMcf(mcf2, mcf2.run(pr), gr, l1, u, c, s2, true, 7620, "#G2");
2.292 - checkMcf(mcf3, mcf3.run(pr), gr, l1, u, c, s3, true, 0, "#G3");
2.293 - checkMcf(mcf4, mcf4.run(pr), gr, l2, u, c, s1, true, 5970, "#G4");
2.294 - checkMcf(mcf5, mcf5.run(pr), gr, l2, u, c, s2, true, 8010, "#G5");
2.295 - checkMcf(mcf6, mcf6.run(pr), gr, l2, u, c, s3, false, 0, "#G6");
2.296 - }
2.297 -
2.298 - // H. Test NetworkSimplex with CANDIDATE_LIST_PIVOT
2.299 - {
2.300 - NetworkSimplex<Digraph>::PivotRuleEnum pr =
2.301 - NetworkSimplex<Digraph>::CANDIDATE_LIST_PIVOT;
2.302 - NetworkSimplex<Digraph> mcf1(gr, u, c, s1);
2.303 - NetworkSimplex<Digraph> mcf2(gr, u, c, v, w, 27);
2.304 - NetworkSimplex<Digraph> mcf3(gr, u, c, s3);
2.305 - NetworkSimplex<Digraph> mcf4(gr, l2, u, c, s1);
2.306 - NetworkSimplex<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.307 - NetworkSimplex<Digraph> mcf6(gr, l2, u, c, s3);
2.308 -
2.309 - checkMcf(mcf1, mcf1.run(pr), gr, l1, u, c, s1, true, 5240, "#H1");
2.310 - checkMcf(mcf2, mcf2.run(pr), gr, l1, u, c, s2, true, 7620, "#H2");
2.311 - checkMcf(mcf3, mcf3.run(pr), gr, l1, u, c, s3, true, 0, "#H3");
2.312 - checkMcf(mcf4, mcf4.run(pr), gr, l2, u, c, s1, true, 5970, "#H4");
2.313 - checkMcf(mcf5, mcf5.run(pr), gr, l2, u, c, s2, true, 8010, "#H5");
2.314 - checkMcf(mcf6, mcf6.run(pr), gr, l2, u, c, s3, false, 0, "#H6");
2.315 - }
2.316 -
2.317 - // I. Test NetworkSimplex with ALTERING_LIST_PIVOT
2.318 - {
2.319 - NetworkSimplex<Digraph>::PivotRuleEnum pr =
2.320 - NetworkSimplex<Digraph>::ALTERING_LIST_PIVOT;
2.321 - NetworkSimplex<Digraph> mcf1(gr, u, c, s1);
2.322 - NetworkSimplex<Digraph> mcf2(gr, u, c, v, w, 27);
2.323 - NetworkSimplex<Digraph> mcf3(gr, u, c, s3);
2.324 - NetworkSimplex<Digraph> mcf4(gr, l2, u, c, s1);
2.325 - NetworkSimplex<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.326 - NetworkSimplex<Digraph> mcf6(gr, l2, u, c, s3);
2.327 -
2.328 - checkMcf(mcf1, mcf1.run(pr), gr, l1, u, c, s1, true, 5240, "#I1");
2.329 - checkMcf(mcf2, mcf2.run(pr), gr, l1, u, c, s2, true, 7620, "#I2");
2.330 - checkMcf(mcf3, mcf3.run(pr), gr, l1, u, c, s3, true, 0, "#I3");
2.331 - checkMcf(mcf4, mcf4.run(pr), gr, l2, u, c, s1, true, 5970, "#I4");
2.332 - checkMcf(mcf5, mcf5.run(pr), gr, l2, u, c, s2, true, 8010, "#I5");
2.333 - checkMcf(mcf6, mcf6.run(pr), gr, l2, u, c, s3, false, 0, "#I6");
2.334 - }
2.335 -
2.336 -/*
2.337 - // J. Test MinCostFlow
2.338 - {
2.339 - MinCostFlow<Digraph> mcf1(gr, u, c, s1);
2.340 - MinCostFlow<Digraph> mcf2(gr, u, c, v, w, 27);
2.341 - MinCostFlow<Digraph> mcf3(gr, u, c, s3);
2.342 - MinCostFlow<Digraph> mcf4(gr, l2, u, c, s1);
2.343 - MinCostFlow<Digraph> mcf5(gr, l2, u, c, v, w, 27);
2.344 - MinCostFlow<Digraph> mcf6(gr, l2, u, c, s3);
2.345 -
2.346 - checkMcf(mcf1, mcf1.run(), gr, l1, u, c, s1, true, 5240, "#J1");
2.347 - checkMcf(mcf2, mcf2.run(), gr, l1, u, c, s2, true, 7620, "#J2");
2.348 - checkMcf(mcf3, mcf3.run(), gr, l1, u, c, s3, true, 0, "#J3");
2.349 - checkMcf(mcf4, mcf4.run(), gr, l2, u, c, s1, true, 5970, "#J4");
2.350 - checkMcf(mcf5, mcf5.run(), gr, l2, u, c, s2, true, 8010, "#J5");
2.351 - checkMcf(mcf6, mcf6.run(), gr, l2, u, c, s3, false, 0, "#J6");
2.352 - }
2.353 -*/
2.354 -/*
2.355 - // K. Test MinCostMaxFlow
2.356 - {
2.357 - MinCostMaxFlow<Digraph> mcmf(gr, u, c, v, w);
2.358 - mcmf.run();
2.359 - checkMcf(mcmf, true, gr, l1, u, c, s3, true, 7620, "#K1");
2.360 - }
2.361 -*/
2.362 -
2.363 return 0;
2.364 }
3.1 --- a/tools/dimacs-solver.cc Tue Mar 24 00:18:25 2009 +0100
3.2 +++ b/tools/dimacs-solver.cc Wed Mar 25 15:58:44 2009 +0100
3.3 @@ -105,9 +105,8 @@
3.4 readDimacsMin(is, g, lower, cap, cost, sup, desc);
3.5 if (report) std::cerr << "Read the file: " << ti << '\n';
3.6 ti.restart();
3.7 - NetworkSimplex< Digraph, Digraph::ArcMap<Value>, Digraph::ArcMap<Value>,
3.8 - Digraph::ArcMap<Value>, Digraph::NodeMap<Value> >
3.9 - ns(g, lower, cap, cost, sup);
3.10 + NetworkSimplex<Digraph, Value> ns(g);
3.11 + ns.lowerMap(lower).capacityMap(cap).costMap(cost).supplyMap(sup);
3.12 if (report) std::cerr << "Setup NetworkSimplex class: " << ti << '\n';
3.13 ti.restart();
3.14 ns.run();