deba@2440: /* -*- C++ -*- deba@2440: * deba@2440: * This file is a part of LEMON, a generic C++ optimization library deba@2440: * deba@2440: * Copyright (C) 2003-2007 deba@2440: * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport deba@2440: * (Egervary Research Group on Combinatorial Optimization, EGRES). deba@2440: * deba@2440: * Permission to use, modify and distribute this software is granted deba@2440: * provided that this copyright notice appears in all copies. For deba@2440: * precise terms see the accompanying LICENSE file. deba@2440: * deba@2440: * This software is provided "AS IS" with no warranty of any kind, deba@2440: * express or implied, and with no claim as to its suitability for any deba@2440: * purpose. deba@2440: * deba@2440: */ deba@2440: deba@2440: #ifndef LEMON_CAPACITY_SCALING_H deba@2440: #define LEMON_CAPACITY_SCALING_H deba@2440: deba@2440: /// \ingroup min_cost_flow deba@2440: /// deba@2440: /// \file deba@2440: /// \brief The capacity scaling algorithm for finding a minimum cost deba@2440: /// flow. deba@2440: kpeter@2535: #include kpeter@2535: #include deba@2440: #include deba@2457: deba@2440: namespace lemon { deba@2440: deba@2440: /// \addtogroup min_cost_flow deba@2440: /// @{ deba@2440: deba@2440: /// \brief Implementation of the capacity scaling version of the kpeter@2509: /// successive shortest path algorithm for finding a minimum cost kpeter@2509: /// flow. deba@2440: /// kpeter@2535: /// \ref CapacityScaling implements the capacity scaling version kpeter@2535: /// of the successive shortest path algorithm for finding a minimum kpeter@2535: /// cost flow. deba@2440: /// deba@2440: /// \param Graph The directed graph type the algorithm runs on. deba@2440: /// \param LowerMap The type of the lower bound map. deba@2440: /// \param CapacityMap The type of the capacity (upper bound) map. deba@2440: /// \param CostMap The type of the cost (length) map. deba@2440: /// \param SupplyMap The type of the supply map. deba@2440: /// deba@2440: /// \warning deba@2440: /// - Edge capacities and costs should be nonnegative integers. kpeter@2535: /// However \c CostMap::Value should be signed type. kpeter@2509: /// - Supply values should be signed integers. deba@2440: /// - \c LowerMap::Value must be convertible to kpeter@2535: /// \c CapacityMap::Value and \c CapacityMap::Value must be kpeter@2535: /// convertible to \c SupplyMap::Value. deba@2440: /// deba@2440: /// \author Peter Kovacs deba@2440: kpeter@2533: template < typename Graph, kpeter@2535: typename LowerMap = typename Graph::template EdgeMap, kpeter@2535: typename CapacityMap = LowerMap, kpeter@2535: typename CostMap = typename Graph::template EdgeMap, kpeter@2535: typename SupplyMap = typename Graph::template NodeMap kpeter@2535: > deba@2440: class CapacityScaling deba@2440: { deba@2440: typedef typename Graph::Node Node; deba@2440: typedef typename Graph::NodeIt NodeIt; deba@2440: typedef typename Graph::Edge Edge; deba@2440: typedef typename Graph::EdgeIt EdgeIt; deba@2440: typedef typename Graph::InEdgeIt InEdgeIt; deba@2440: typedef typename Graph::OutEdgeIt OutEdgeIt; deba@2440: deba@2440: typedef typename LowerMap::Value Lower; deba@2440: typedef typename CapacityMap::Value Capacity; deba@2440: typedef typename CostMap::Value Cost; deba@2440: typedef typename SupplyMap::Value Supply; deba@2440: typedef typename Graph::template EdgeMap CapacityRefMap; deba@2440: typedef typename Graph::template NodeMap SupplyRefMap; kpeter@2535: typedef typename Graph::template NodeMap PredMap; deba@2440: deba@2440: public: deba@2440: kpeter@2535: /// \brief Type to enable or disable capacity scaling. kpeter@2535: enum ScalingEnum { kpeter@2535: WITH_SCALING = 0, kpeter@2535: WITHOUT_SCALING = -1 kpeter@2535: }; kpeter@2535: deba@2440: /// \brief The type of the flow map. deba@2440: typedef CapacityRefMap FlowMap; deba@2440: /// \brief The type of the potential map. deba@2440: typedef typename Graph::template NodeMap PotentialMap; deba@2440: deba@2440: protected: deba@2440: kpeter@2535: /// \brief Special implementation of the \ref Dijkstra algorithm kpeter@2535: /// for finding shortest paths in the residual network of the graph kpeter@2535: /// with respect to the reduced edge costs and modifying the kpeter@2535: /// node potentials according to the distance of the nodes. kpeter@2535: class ResidualDijkstra deba@2440: { kpeter@2535: typedef typename Graph::template NodeMap CostNodeMap; kpeter@2535: typedef typename Graph::template NodeMap PredMap; deba@2440: kpeter@2535: typedef typename Graph::template NodeMap HeapCrossRef; kpeter@2535: typedef BinHeap Heap; kpeter@2535: kpeter@2535: protected: kpeter@2535: kpeter@2535: /// \brief The directed graph the algorithm runs on. kpeter@2535: const Graph &graph; kpeter@2535: kpeter@2535: /// \brief The flow map. kpeter@2535: const FlowMap &flow; kpeter@2535: /// \brief The residual capacity map. kpeter@2535: const CapacityRefMap &res_cap; kpeter@2535: /// \brief The cost map. kpeter@2535: const CostMap &cost; kpeter@2535: /// \brief The excess map. kpeter@2535: const SupplyRefMap &excess; kpeter@2535: kpeter@2535: /// \brief The potential map. kpeter@2535: PotentialMap &potential; kpeter@2535: kpeter@2535: /// \brief The distance map. kpeter@2535: CostNodeMap dist; kpeter@2535: /// \brief The map of predecessors edges. kpeter@2535: PredMap &pred; kpeter@2535: /// \brief The processed (i.e. permanently labeled) nodes. kpeter@2535: std::vector proc_nodes; deba@2440: deba@2440: public: deba@2440: kpeter@2535: /// \brief The constructor of the class. kpeter@2535: ResidualDijkstra( const Graph &_graph, kpeter@2535: const FlowMap &_flow, kpeter@2535: const CapacityRefMap &_res_cap, kpeter@2535: const CostMap &_cost, kpeter@2535: const SupplyMap &_excess, kpeter@2535: PotentialMap &_potential, kpeter@2535: PredMap &_pred ) : kpeter@2535: graph(_graph), flow(_flow), res_cap(_res_cap), cost(_cost), kpeter@2535: excess(_excess), potential(_potential), dist(_graph), kpeter@2535: pred(_pred) kpeter@2535: {} deba@2440: kpeter@2535: /// \brief Runs the algorithm from the given source node. kpeter@2535: Node run(Node s, Capacity delta) { kpeter@2535: HeapCrossRef heap_cross_ref(graph, Heap::PRE_HEAP); kpeter@2535: Heap heap(heap_cross_ref); kpeter@2535: heap.push(s, 0); kpeter@2535: pred[s] = INVALID; kpeter@2535: proc_nodes.clear(); kpeter@2535: kpeter@2535: // Processing nodes kpeter@2535: while (!heap.empty() && excess[heap.top()] > -delta) { kpeter@2535: Node u = heap.top(), v; kpeter@2535: Cost d = heap.prio() - potential[u], nd; kpeter@2535: dist[u] = heap.prio(); kpeter@2535: heap.pop(); kpeter@2535: proc_nodes.push_back(u); kpeter@2535: kpeter@2535: // Traversing outgoing edges kpeter@2535: for (OutEdgeIt e(graph, u); e != INVALID; ++e) { kpeter@2535: if (res_cap[e] >= delta) { kpeter@2535: v = graph.target(e); kpeter@2535: switch(heap.state(v)) { kpeter@2535: case Heap::PRE_HEAP: kpeter@2535: heap.push(v, d + cost[e] + potential[v]); kpeter@2535: pred[v] = e; kpeter@2535: break; kpeter@2535: case Heap::IN_HEAP: kpeter@2535: nd = d + cost[e] + potential[v]; kpeter@2535: if (nd < heap[v]) { kpeter@2535: heap.decrease(v, nd); kpeter@2535: pred[v] = e; kpeter@2535: } kpeter@2535: break; kpeter@2535: case Heap::POST_HEAP: kpeter@2535: break; kpeter@2535: } kpeter@2535: } kpeter@2535: } kpeter@2535: kpeter@2535: // Traversing incoming edges kpeter@2535: for (InEdgeIt e(graph, u); e != INVALID; ++e) { kpeter@2535: if (flow[e] >= delta) { kpeter@2535: v = graph.source(e); kpeter@2535: switch(heap.state(v)) { kpeter@2535: case Heap::PRE_HEAP: kpeter@2535: heap.push(v, d - cost[e] + potential[v]); kpeter@2535: pred[v] = e; kpeter@2535: break; kpeter@2535: case Heap::IN_HEAP: kpeter@2535: nd = d - cost[e] + potential[v]; kpeter@2535: if (nd < heap[v]) { kpeter@2535: heap.decrease(v, nd); kpeter@2535: pred[v] = e; kpeter@2535: } kpeter@2535: break; kpeter@2535: case Heap::POST_HEAP: kpeter@2535: break; kpeter@2535: } kpeter@2535: } kpeter@2535: } kpeter@2535: } kpeter@2535: if (heap.empty()) return INVALID; kpeter@2535: kpeter@2535: // Updating potentials of processed nodes kpeter@2535: Node t = heap.top(); kpeter@2535: Cost dt = heap.prio(); kpeter@2535: for (int i = 0; i < proc_nodes.size(); ++i) kpeter@2535: potential[proc_nodes[i]] -= dist[proc_nodes[i]] - dt; kpeter@2535: kpeter@2535: return t; deba@2440: } deba@2440: kpeter@2535: }; //class ResidualDijkstra deba@2440: deba@2440: protected: deba@2440: deba@2440: /// \brief The directed graph the algorithm runs on. deba@2440: const Graph &graph; deba@2440: /// \brief The original lower bound map. deba@2440: const LowerMap *lower; deba@2440: /// \brief The modified capacity map. deba@2440: CapacityRefMap capacity; deba@2440: /// \brief The cost map. deba@2440: const CostMap &cost; deba@2440: /// \brief The modified supply map. deba@2440: SupplyRefMap supply; deba@2440: /// \brief The sum of supply values equals zero. deba@2440: bool valid_supply; deba@2440: deba@2440: /// \brief The edge map of the current flow. deba@2440: FlowMap flow; deba@2440: /// \brief The potential node map. deba@2440: PotentialMap potential; deba@2440: kpeter@2535: /// \brief The residual capacity map. kpeter@2535: CapacityRefMap res_cap; kpeter@2535: /// \brief The excess map. kpeter@2535: SupplyRefMap excess; kpeter@2535: /// \brief The excess nodes (i.e. nodes with positive excess). kpeter@2535: std::vector excess_nodes; deba@2440: /// \brief The index of the next excess node. deba@2440: int next_node; deba@2440: kpeter@2535: /// \brief The scaling status (enabled or disabled). kpeter@2535: ScalingEnum scaling; deba@2440: /// \brief The delta parameter used for capacity scaling. deba@2440: Capacity delta; kpeter@2535: /// \brief The maximum number of phases. kpeter@2535: Capacity phase_num; kpeter@2535: /// \brief The deficit nodes. kpeter@2535: std::vector deficit_nodes; deba@2440: kpeter@2535: /// \brief Implementation of the \ref Dijkstra algorithm for kpeter@2535: /// finding augmenting shortest paths in the residual network. kpeter@2535: ResidualDijkstra dijkstra; kpeter@2535: /// \brief The map of predecessors edges. kpeter@2535: PredMap pred; deba@2440: deba@2440: public : deba@2440: deba@2440: /// \brief General constructor of the class (with lower bounds). deba@2440: /// deba@2440: /// General constructor of the class (with lower bounds). deba@2440: /// deba@2440: /// \param _graph The directed graph the algorithm runs on. deba@2440: /// \param _lower The lower bounds of the edges. deba@2440: /// \param _capacity The capacities (upper bounds) of the edges. deba@2440: /// \param _cost The cost (length) values of the edges. deba@2440: /// \param _supply The supply values of the nodes (signed). deba@2440: CapacityScaling( const Graph &_graph, kpeter@2535: const LowerMap &_lower, kpeter@2535: const CapacityMap &_capacity, kpeter@2535: const CostMap &_cost, kpeter@2535: const SupplyMap &_supply ) : deba@2440: graph(_graph), lower(&_lower), capacity(_graph), cost(_cost), deba@2440: supply(_graph), flow(_graph, 0), potential(_graph, 0), kpeter@2535: res_cap(_graph), excess(_graph), pred(_graph), kpeter@2535: dijkstra(graph, flow, res_cap, cost, excess, potential, pred) deba@2440: { deba@2440: // Removing nonzero lower bounds deba@2440: capacity = subMap(_capacity, _lower); kpeter@2535: res_cap = capacity; deba@2440: Supply sum = 0; deba@2440: for (NodeIt n(graph); n != INVALID; ++n) { kpeter@2535: Supply s = _supply[n]; kpeter@2535: for (InEdgeIt e(graph, n); e != INVALID; ++e) kpeter@2535: s += _lower[e]; kpeter@2535: for (OutEdgeIt e(graph, n); e != INVALID; ++e) kpeter@2535: s -= _lower[e]; kpeter@2535: supply[n] = s; kpeter@2535: sum += s; deba@2440: } deba@2440: valid_supply = sum == 0; deba@2440: } deba@2440: deba@2440: /// \brief General constructor of the class (without lower bounds). deba@2440: /// deba@2440: /// General constructor of the class (without lower bounds). deba@2440: /// deba@2440: /// \param _graph The directed graph the algorithm runs on. deba@2440: /// \param _capacity The capacities (upper bounds) of the edges. deba@2440: /// \param _cost The cost (length) values of the edges. deba@2440: /// \param _supply The supply values of the nodes (signed). deba@2440: CapacityScaling( const Graph &_graph, kpeter@2535: const CapacityMap &_capacity, kpeter@2535: const CostMap &_cost, kpeter@2535: const SupplyMap &_supply ) : deba@2440: graph(_graph), lower(NULL), capacity(_capacity), cost(_cost), deba@2440: supply(_supply), flow(_graph, 0), potential(_graph, 0), kpeter@2535: res_cap(_capacity), excess(_graph), pred(_graph), kpeter@2535: dijkstra(graph, flow, res_cap, cost, excess, potential) deba@2440: { deba@2440: // Checking the sum of supply values deba@2440: Supply sum = 0; deba@2440: for (NodeIt n(graph); n != INVALID; ++n) sum += supply[n]; deba@2440: valid_supply = sum == 0; deba@2440: } deba@2440: deba@2440: /// \brief Simple constructor of the class (with lower bounds). deba@2440: /// deba@2440: /// Simple constructor of the class (with lower bounds). deba@2440: /// deba@2440: /// \param _graph The directed graph the algorithm runs on. deba@2440: /// \param _lower The lower bounds of the edges. deba@2440: /// \param _capacity The capacities (upper bounds) of the edges. deba@2440: /// \param _cost The cost (length) values of the edges. deba@2440: /// \param _s The source node. deba@2440: /// \param _t The target node. deba@2440: /// \param _flow_value The required amount of flow from node \c _s deba@2440: /// to node \c _t (i.e. the supply of \c _s and the demand of deba@2440: /// \c _t). deba@2440: CapacityScaling( const Graph &_graph, kpeter@2535: const LowerMap &_lower, kpeter@2535: const CapacityMap &_capacity, kpeter@2535: const CostMap &_cost, kpeter@2535: Node _s, Node _t, kpeter@2535: Supply _flow_value ) : deba@2440: graph(_graph), lower(&_lower), capacity(_graph), cost(_cost), deba@2440: supply(_graph), flow(_graph, 0), potential(_graph, 0), kpeter@2535: res_cap(_graph), excess(_graph), pred(_graph), kpeter@2535: dijkstra(graph, flow, res_cap, cost, excess, potential) deba@2440: { deba@2440: // Removing nonzero lower bounds deba@2440: capacity = subMap(_capacity, _lower); kpeter@2535: res_cap = capacity; deba@2440: for (NodeIt n(graph); n != INVALID; ++n) { kpeter@2535: Supply s = 0; kpeter@2535: if (n == _s) s = _flow_value; kpeter@2535: if (n == _t) s = -_flow_value; kpeter@2535: for (InEdgeIt e(graph, n); e != INVALID; ++e) kpeter@2535: s += _lower[e]; kpeter@2535: for (OutEdgeIt e(graph, n); e != INVALID; ++e) kpeter@2535: s -= _lower[e]; kpeter@2535: supply[n] = s; deba@2440: } deba@2440: valid_supply = true; deba@2440: } deba@2440: deba@2440: /// \brief Simple constructor of the class (without lower bounds). deba@2440: /// deba@2440: /// Simple constructor of the class (without lower bounds). deba@2440: /// deba@2440: /// \param _graph The directed graph the algorithm runs on. deba@2440: /// \param _capacity The capacities (upper bounds) of the edges. deba@2440: /// \param _cost The cost (length) values of the edges. deba@2440: /// \param _s The source node. deba@2440: /// \param _t The target node. deba@2440: /// \param _flow_value The required amount of flow from node \c _s deba@2440: /// to node \c _t (i.e. the supply of \c _s and the demand of deba@2440: /// \c _t). deba@2440: CapacityScaling( const Graph &_graph, kpeter@2535: const CapacityMap &_capacity, kpeter@2535: const CostMap &_cost, kpeter@2535: Node _s, Node _t, kpeter@2535: Supply _flow_value ) : deba@2440: graph(_graph), lower(NULL), capacity(_capacity), cost(_cost), deba@2440: supply(_graph, 0), flow(_graph, 0), potential(_graph, 0), kpeter@2535: res_cap(_capacity), excess(_graph), pred(_graph), kpeter@2535: dijkstra(graph, flow, res_cap, cost, excess, potential) deba@2440: { deba@2440: supply[_s] = _flow_value; deba@2440: supply[_t] = -_flow_value; deba@2440: valid_supply = true; deba@2440: } deba@2440: deba@2440: /// \brief Returns a const reference to the flow map. deba@2440: /// deba@2440: /// Returns a const reference to the flow map. deba@2440: /// deba@2440: /// \pre \ref run() must be called before using this function. deba@2440: const FlowMap& flowMap() const { deba@2440: return flow; deba@2440: } deba@2440: deba@2440: /// \brief Returns a const reference to the potential map (the dual deba@2440: /// solution). deba@2440: /// deba@2440: /// Returns a const reference to the potential map (the dual deba@2440: /// solution). deba@2440: /// deba@2440: /// \pre \ref run() must be called before using this function. deba@2440: const PotentialMap& potentialMap() const { deba@2440: return potential; deba@2440: } deba@2440: deba@2440: /// \brief Returns the total cost of the found flow. deba@2440: /// deba@2440: /// Returns the total cost of the found flow. The complexity of the deba@2440: /// function is \f$ O(e) \f$. deba@2440: /// deba@2440: /// \pre \ref run() must be called before using this function. deba@2440: Cost totalCost() const { deba@2440: Cost c = 0; deba@2440: for (EdgeIt e(graph); e != INVALID; ++e) kpeter@2535: c += flow[e] * cost[e]; deba@2440: return c; deba@2440: } deba@2440: deba@2457: /// \brief Runs the algorithm. deba@2440: /// deba@2457: /// Runs the algorithm. deba@2440: /// kpeter@2535: /// \param scaling_mode The scaling mode. In case of WITH_SCALING kpeter@2535: /// capacity scaling is enabled in the algorithm (this is the kpeter@2535: /// default value) otherwise it is disabled. kpeter@2535: /// If the maximum edge capacity and/or the amount of total supply kpeter@2535: /// is small, the algorithm could be faster without scaling. kpeter@2535: /// deba@2440: /// \return \c true if a feasible flow can be found. kpeter@2535: bool run(int scaling_mode = WITH_SCALING) { kpeter@2535: return init(scaling_mode) && start(); deba@2440: } deba@2440: deba@2440: protected: deba@2440: deba@2440: /// \brief Initializes the algorithm. kpeter@2535: bool init(int scaling_mode) { deba@2440: if (!valid_supply) return false; kpeter@2535: excess = supply; deba@2440: deba@2440: // Initilaizing delta value kpeter@2535: if (scaling_mode == WITH_SCALING) { kpeter@2535: // With scaling kpeter@2535: Supply max_sup = 0, max_dem = 0; kpeter@2535: for (NodeIt n(graph); n != INVALID; ++n) { kpeter@2535: if ( supply[n] > max_sup) max_sup = supply[n]; kpeter@2535: if (-supply[n] > max_dem) max_dem = -supply[n]; kpeter@2535: } kpeter@2535: if (max_dem < max_sup) max_sup = max_dem; kpeter@2535: phase_num = 0; kpeter@2535: for (delta = 1; 2 * delta <= max_sup; delta *= 2) kpeter@2535: ++phase_num; kpeter@2535: } else { kpeter@2535: // Without scaling kpeter@2535: delta = 1; deba@2440: } deba@2440: return true; deba@2440: } deba@2440: kpeter@2535: /// \brief Executes the algorithm. kpeter@2535: bool start() { kpeter@2535: if (delta > 1) kpeter@2535: return startWithScaling(); kpeter@2535: else kpeter@2535: return startWithoutScaling(); kpeter@2535: } kpeter@2535: kpeter@2535: /// \brief Executes the capacity scaling version of the successive kpeter@2535: /// shortest path algorithm. kpeter@2535: bool startWithScaling() { kpeter@2535: // Processing capacity scaling phases kpeter@2535: Node s, t; kpeter@2535: int phase_cnt = 0; kpeter@2535: int factor = 4; kpeter@2535: while (true) { kpeter@2535: // Saturating all edges not satisfying the optimality condition kpeter@2535: for (EdgeIt e(graph); e != INVALID; ++e) { kpeter@2535: Node u = graph.source(e), v = graph.target(e); kpeter@2535: Cost c = cost[e] - potential[u] + potential[v]; kpeter@2535: if (c < 0 && res_cap[e] >= delta) { kpeter@2535: excess[u] -= res_cap[e]; kpeter@2535: excess[v] += res_cap[e]; kpeter@2535: flow[e] = capacity[e]; kpeter@2535: res_cap[e] = 0; kpeter@2535: } kpeter@2535: else if (c > 0 && flow[e] >= delta) { kpeter@2535: excess[u] += flow[e]; kpeter@2535: excess[v] -= flow[e]; kpeter@2535: flow[e] = 0; kpeter@2535: res_cap[e] = capacity[e]; kpeter@2535: } kpeter@2535: } kpeter@2535: kpeter@2535: // Finding excess nodes and deficit nodes kpeter@2535: excess_nodes.clear(); kpeter@2535: deficit_nodes.clear(); kpeter@2535: for (NodeIt n(graph); n != INVALID; ++n) { kpeter@2535: if (excess[n] >= delta) excess_nodes.push_back(n); kpeter@2535: if (excess[n] <= -delta) deficit_nodes.push_back(n); kpeter@2535: } kpeter@2535: next_node = 0; kpeter@2535: kpeter@2535: // Finding augmenting shortest paths kpeter@2535: while (next_node < excess_nodes.size()) { kpeter@2535: // Checking deficit nodes kpeter@2535: if (delta > 1) { kpeter@2535: bool delta_deficit = false; kpeter@2535: for (int i = 0; i < deficit_nodes.size(); ++i) { kpeter@2535: if (excess[deficit_nodes[i]] <= -delta) { kpeter@2535: delta_deficit = true; kpeter@2535: break; kpeter@2535: } kpeter@2535: } kpeter@2535: if (!delta_deficit) break; kpeter@2535: } kpeter@2535: kpeter@2535: // Running Dijkstra kpeter@2535: s = excess_nodes[next_node]; kpeter@2535: if ((t = dijkstra.run(s, delta)) == INVALID) { kpeter@2535: if (delta > 1) { kpeter@2535: ++next_node; kpeter@2535: continue; kpeter@2535: } kpeter@2535: return false; kpeter@2535: } kpeter@2535: kpeter@2535: // Augmenting along a shortest path from s to t. kpeter@2535: Capacity d = excess[s] < -excess[t] ? excess[s] : -excess[t]; kpeter@2535: Node u = t; kpeter@2535: Edge e; kpeter@2535: if (d > delta) { kpeter@2535: while ((e = pred[u]) != INVALID) { kpeter@2535: Capacity rc; kpeter@2535: if (u == graph.target(e)) { kpeter@2535: rc = res_cap[e]; kpeter@2535: u = graph.source(e); kpeter@2535: } else { kpeter@2535: rc = flow[e]; kpeter@2535: u = graph.target(e); kpeter@2535: } kpeter@2535: if (rc < d) d = rc; kpeter@2535: } kpeter@2535: } kpeter@2535: u = t; kpeter@2535: while ((e = pred[u]) != INVALID) { kpeter@2535: if (u == graph.target(e)) { kpeter@2535: flow[e] += d; kpeter@2535: res_cap[e] -= d; kpeter@2535: u = graph.source(e); kpeter@2535: } else { kpeter@2535: flow[e] -= d; kpeter@2535: res_cap[e] += d; kpeter@2535: u = graph.target(e); kpeter@2535: } kpeter@2535: } kpeter@2535: excess[s] -= d; kpeter@2535: excess[t] += d; kpeter@2535: kpeter@2535: if (excess[s] < delta) ++next_node; kpeter@2535: } kpeter@2535: kpeter@2535: if (delta == 1) break; kpeter@2535: if (++phase_cnt > phase_num / 4) factor = 2; kpeter@2535: delta = delta <= factor ? 1 : delta / factor; kpeter@2535: } kpeter@2535: kpeter@2535: // Handling nonzero lower bounds kpeter@2535: if (lower) { kpeter@2535: for (EdgeIt e(graph); e != INVALID; ++e) kpeter@2535: flow[e] += (*lower)[e]; kpeter@2535: } kpeter@2535: return true; kpeter@2535: } kpeter@2535: deba@2440: /// \brief Executes the successive shortest path algorithm without deba@2440: /// capacity scaling. kpeter@2535: bool startWithoutScaling() { deba@2440: // Finding excess nodes kpeter@2535: for (NodeIt n(graph); n != INVALID; ++n) { kpeter@2535: if (excess[n] > 0) excess_nodes.push_back(n); deba@2440: } deba@2440: if (excess_nodes.size() == 0) return true; deba@2440: next_node = 0; deba@2440: deba@2457: // Finding shortest paths kpeter@2535: Node s, t; kpeter@2535: while ( excess[excess_nodes[next_node]] > 0 || kpeter@2535: ++next_node < excess_nodes.size() ) deba@2440: { kpeter@2535: // Running Dijkstra kpeter@2535: s = excess_nodes[next_node]; kpeter@2535: if ((t = dijkstra.run(s, 1)) == INVALID) kpeter@2535: return false; deba@2440: kpeter@2535: // Augmenting along a shortest path from s to t kpeter@2535: Capacity d = excess[s] < -excess[t] ? excess[s] : -excess[t]; kpeter@2535: Node u = t; kpeter@2535: Edge e; kpeter@2535: while ((e = pred[u]) != INVALID) { kpeter@2535: Capacity rc; kpeter@2535: if (u == graph.target(e)) { kpeter@2535: rc = res_cap[e]; kpeter@2535: u = graph.source(e); kpeter@2535: } else { kpeter@2535: rc = flow[e]; kpeter@2535: u = graph.target(e); kpeter@2535: } kpeter@2535: if (rc < d) d = rc; kpeter@2535: } kpeter@2535: u = t; kpeter@2535: while ((e = pred[u]) != INVALID) { kpeter@2535: if (u == graph.target(e)) { kpeter@2535: flow[e] += d; kpeter@2535: res_cap[e] -= d; kpeter@2535: u = graph.source(e); kpeter@2535: } else { kpeter@2535: flow[e] -= d; kpeter@2535: res_cap[e] += d; kpeter@2535: u = graph.target(e); kpeter@2535: } kpeter@2535: } kpeter@2535: excess[s] -= d; kpeter@2535: excess[t] += d; deba@2440: } deba@2440: deba@2440: // Handling nonzero lower bounds deba@2440: if (lower) { kpeter@2535: for (EdgeIt e(graph); e != INVALID; ++e) kpeter@2535: flow[e] += (*lower)[e]; deba@2440: } deba@2440: return true; deba@2440: } deba@2440: deba@2440: }; //class CapacityScaling deba@2440: deba@2440: ///@} deba@2440: deba@2440: } //namespace lemon deba@2440: deba@2440: #endif //LEMON_CAPACITY_SCALING_H