COIN-OR::LEMON - Graph Library

Changeset 2586:37fb2c384c78 in lemon-0.x


Ignore:
Timestamp:
02/29/08 16:55:13 (11 years ago)
Author:
Peter Kovacs
Branch:
default
Phase:
public
Convert:
svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@3468
Message:

Reimplemented Suurballe class.

  • The new version is the specialized version of CapacityScaling?.
  • It is about 10-20 times faster than the former Suurballe algorithm

and about 20-50 percent faster than CapacityScaling?.

  • Doc improvements.
  • The test file is also replaced.
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • lemon/suurballe.h

    r2553 r2586  
    2222///\ingroup shortest_path
    2323///\file
    24 ///\brief An algorithm for finding k paths of minimal total length.
    25 
    26 
    27 #include <lemon/maps.h>
     24///\brief An algorithm for finding edge-disjoint paths between two
     25/// nodes having minimum total length.
     26
    2827#include <vector>
     28#include <lemon/bin_heap.h>
    2929#include <lemon/path.h>
    30 #include <lemon/ssp_min_cost_flow.h>
    3130
    3231namespace lemon {
    3332
    34 /// \addtogroup shortest_path
    35 /// @{
    36 
    37   ///\brief Implementation of an algorithm for finding k edge-disjoint
    38   /// paths between 2 nodes of minimal total length
    39   ///
    40   /// The class \ref lemon::Suurballe implements
    41   /// an algorithm for finding k edge-disjoint paths
    42   /// from a given source node to a given target node in an
    43   /// edge-weighted directed graph having minimal total weight (length).
    44   ///
    45   ///\warning Length values should be nonnegative!
    46   ///
    47   ///\param Graph The directed graph type the algorithm runs on.
    48   ///\param LengthMap The type of the length map (values should be nonnegative).
    49   ///
    50   ///\note It it questionable whether it is correct to call this method after
    51   ///%Suurballe for it is just a special case of Edmonds' and Karp's algorithm
    52   ///for finding minimum cost flows. In fact, this implementation just
    53   ///wraps the SspMinCostFlow algorithms. The paper of both %Suurballe and
    54   ///Edmonds-Karp published in 1972, therefore it is possibly right to
    55   ///state that they are
    56   ///independent results. Most frequently this special case is referred as
    57   ///%Suurballe method in the literature, especially in communication
    58   ///network context.
    59   ///\author Attila Bernath
    60   template <typename Graph, typename LengthMap>
    61   class Suurballe{
    62 
     33  /// \addtogroup shortest_path
     34  /// @{
     35
     36  /// \brief Implementation of an algorithm for finding edge-disjoint
     37  /// paths between two nodes having minimum total length.
     38  ///
     39  /// \ref lemon::Suurballe "Suurballe" implements an algorithm for
     40  /// finding edge-disjoint paths having minimum total length (cost)
     41  /// from a given source node to a given target node in a directed
     42  /// graph.
     43  ///
     44  /// In fact, this implementation is the specialization of the
     45  /// \ref CapacityScaling "successive shortest path" algorithm.
     46  ///
     47  /// \tparam Graph The directed graph type the algorithm runs on.
     48  /// \tparam LengthMap The type of the length (cost) map.
     49  ///
     50  /// \warning Length values should be \e non-negative \e integers.
     51  ///
     52  /// \note For finding node-disjoint paths this algorithm can be used
     53  /// with \ref SplitGraphAdaptor.
     54  ///
     55  /// \author Attila Bernath and Peter Kovacs
     56 
     57  template < typename Graph,
     58             typename LengthMap = typename Graph::template EdgeMap<int> >
     59  class Suurballe
     60  {
     61    GRAPH_TYPEDEFS(typename Graph);
    6362
    6463    typedef typename LengthMap::Value Length;
     64    typedef ConstMap<Edge, int> ConstEdgeMap;
     65    typedef typename Graph::template NodeMap<Edge> PredMap;
     66
     67  public:
     68
     69    /// The type of the flow map.
     70    typedef typename Graph::template EdgeMap<int> FlowMap;
     71    /// The type of the potential map.
     72    typedef typename Graph::template NodeMap<Length> PotentialMap;
     73    /// The type of the path structures.
     74    typedef SimplePath<Graph> Path;
     75
     76  private:
     77 
     78    /// \brief Special implementation of the \ref Dijkstra algorithm
     79    /// for finding shortest paths in the residual network.
     80    ///
     81    /// \ref ResidualDijkstra is a special implementation of the
     82    /// \ref Dijkstra algorithm for finding shortest paths in the
     83    /// residual network of the graph with respect to the reduced edge
     84    /// lengths and modifying the node potentials according to the
     85    /// distance of the nodes.
     86    class ResidualDijkstra
     87    {
     88      typedef typename Graph::template NodeMap<int> HeapCrossRef;
     89      typedef BinHeap<Length, HeapCrossRef> Heap;
     90
     91    private:
     92
     93      // The directed graph the algorithm runs on
     94      const Graph &_graph;
     95
     96      // The main maps
     97      const FlowMap &_flow;
     98      const LengthMap &_length;
     99      PotentialMap &_potential;
     100
     101      // The distance map
     102      PotentialMap _dist;
     103      // The pred edge map
     104      PredMap &_pred;
     105      // The processed (i.e. permanently labeled) nodes
     106      std::vector<Node> _proc_nodes;
     107     
     108      Node _s;
     109      Node _t;
     110
     111    public:
     112
     113      /// Constructor.
     114      ResidualDijkstra( const Graph &graph,
     115                        const FlowMap &flow,
     116                        const LengthMap &length,
     117                        PotentialMap &potential,
     118                        PredMap &pred,
     119                        Node s, Node t ) :
     120        _graph(graph), _flow(flow), _length(length), _potential(potential),
     121        _dist(graph), _pred(pred), _s(s), _t(t) {}
     122
     123      /// \brief Runs the algorithm. Returns \c true if a path is found
     124      /// from the source node to the target node.
     125      bool run() {
     126        HeapCrossRef heap_cross_ref(_graph, Heap::PRE_HEAP);
     127        Heap heap(heap_cross_ref);
     128        heap.push(_s, 0);
     129        _pred[_s] = INVALID;
     130        _proc_nodes.clear();
     131
     132        // Processing nodes
     133        while (!heap.empty() && heap.top() != _t) {
     134          Node u = heap.top(), v;
     135          Length d = heap.prio() + _potential[u], nd;
     136          _dist[u] = heap.prio();
     137          heap.pop();
     138          _proc_nodes.push_back(u);
     139
     140          // Traversing outgoing edges
     141          for (OutEdgeIt e(_graph, u); e != INVALID; ++e) {
     142            if (_flow[e] == 0) {
     143              v = _graph.target(e);
     144              switch(heap.state(v)) {
     145              case Heap::PRE_HEAP:
     146                heap.push(v, d + _length[e] - _potential[v]);
     147                _pred[v] = e;
     148                break;
     149              case Heap::IN_HEAP:
     150                nd = d + _length[e] - _potential[v];
     151                if (nd < heap[v]) {
     152                  heap.decrease(v, nd);
     153                  _pred[v] = e;
     154                }
     155                break;
     156              case Heap::POST_HEAP:
     157                break;
     158              }
     159            }
     160          }
     161
     162          // Traversing incoming edges
     163          for (InEdgeIt e(_graph, u); e != INVALID; ++e) {
     164            if (_flow[e] == 1) {
     165              v = _graph.source(e);
     166              switch(heap.state(v)) {
     167              case Heap::PRE_HEAP:
     168                heap.push(v, d - _length[e] - _potential[v]);
     169                _pred[v] = e;
     170                break;
     171              case Heap::IN_HEAP:
     172                nd = d - _length[e] - _potential[v];
     173                if (nd < heap[v]) {
     174                  heap.decrease(v, nd);
     175                  _pred[v] = e;
     176                }
     177                break;
     178              case Heap::POST_HEAP:
     179                break;
     180              }
     181            }
     182          }
     183        }
     184        if (heap.empty()) return false;
     185
     186        // Updating potentials of processed nodes
     187        Length t_dist = heap.prio();
     188        for (int i = 0; i < int(_proc_nodes.size()); ++i)
     189          _potential[_proc_nodes[i]] += _dist[_proc_nodes[i]] - t_dist;
     190        return true;
     191      }
     192
     193    }; //class ResidualDijkstra
     194
     195  private:
     196
     197    // The directed graph the algorithm runs on
     198    const Graph &_graph;
     199    // The length map
     200    const LengthMap &_length;
    65201   
    66     typedef typename Graph::Node Node;
    67     typedef typename Graph::NodeIt NodeIt;
    68     typedef typename Graph::Edge Edge;
    69     typedef typename Graph::OutEdgeIt OutEdgeIt;
    70     typedef typename Graph::template EdgeMap<int> EdgeIntMap;
    71 
    72     typedef ConstMap<Edge,int> ConstMap;
    73 
    74     const Graph& G;
    75 
    76     Node s;
    77     Node t;
    78 
    79     //Auxiliary variables
    80     //This is the capacity map for the mincostflow problem
    81     ConstMap const1map;
    82     //This MinCostFlow instance will actually solve the problem
    83     SspMinCostFlow<Graph, LengthMap, ConstMap> min_cost_flow;
    84 
    85     //Container to store found paths
    86     std::vector<SimplePath<Graph> > paths;
    87 
    88   public :
    89 
    90 
    91     /// \brief The constructor of the class.
    92     ///
    93     /// \param _G The directed graph the algorithm runs on.
    94     /// \param _length The length (weight or cost) of the edges.
    95     /// \param _s Source node.
    96     /// \param _t Target node.
    97     Suurballe(Graph& _G, LengthMap& _length, Node _s, Node _t) :
    98       G(_G), s(_s), t(_t), const1map(1),
    99       min_cost_flow(_G, _length, const1map, _s, _t) { }
     202    // Edge map of the current flow
     203    FlowMap *_flow;
     204    bool _local_flow;
     205    // Node map of the current potentials
     206    PotentialMap *_potential;
     207    bool _local_potential;
     208
     209    // The source node
     210    Node _source;
     211    // The target node
     212    Node _target;
     213
     214    // Container to store the found paths
     215    std::vector< SimplePath<Graph> > paths;
     216    int _path_num;
     217
     218    // The pred edge map
     219    PredMap _pred;
     220    // Implementation of the Dijkstra algorithm for finding augmenting
     221    // shortest paths in the residual network
     222    ResidualDijkstra *_dijkstra;
     223
     224  public:
     225
     226    /// \brief Constructor.
     227    ///
     228    /// Constructor.
     229    ///
     230    /// \param graph The directed graph the algorithm runs on.
     231    /// \param length The length (cost) values of the edges.
     232    /// \param s The source node.
     233    /// \param t The target node.
     234    Suurballe( const Graph &graph,
     235               const LengthMap &length,
     236               Node s, Node t ) :
     237      _graph(graph), _length(length), _flow(0), _local_flow(false),
     238      _potential(0), _local_potential(false), _source(s), _target(t),
     239      _pred(graph) {}
     240
     241    /// Destructor.
     242    ~Suurballe() {
     243      if (_local_flow) delete _flow;
     244      if (_local_potential) delete _potential;
     245      delete _dijkstra;
     246    }
     247
     248    /// \brief Sets the flow map.
     249    ///
     250    /// Sets the flow map.
     251    ///
     252    /// The found flow contains only 0 and 1 values. It is the union of
     253    /// the found edge-disjoint paths.
     254    ///
     255    /// \return \c (*this)
     256    Suurballe& flowMap(FlowMap &map) {
     257      if (_local_flow) {
     258        delete _flow;
     259        _local_flow = false;
     260      }
     261      _flow = &map;
     262      return *this;
     263    }
     264
     265    /// \brief Sets the potential map.
     266    ///
     267    /// Sets the potential map.
     268    ///
     269    /// The potentials provide the dual solution of the underlying
     270    /// minimum cost flow problem.
     271    ///
     272    /// \return \c (*this)
     273    Suurballe& potentialMap(PotentialMap &map) {
     274      if (_local_potential) {
     275        delete _potential;
     276        _local_potential = false;
     277      }
     278      _potential = &map;
     279      return *this;
     280    }
     281
     282    /// \name Execution control
     283    /// The simplest way to execute the algorithm is to call the run()
     284    /// function.
     285    /// \n
     286    /// If you only need the flow that is the union of the found
     287    /// edge-disjoint paths, you may call init() and findFlow().
     288
     289    /// @{
    100290
    101291    /// \brief Runs the algorithm.
    102292    ///
    103293    /// Runs the algorithm.
    104     /// Returns k if there are at least k edge-disjoint paths from s to t.
    105     /// Otherwise it returns the number of edge-disjoint paths found
    106     /// from s to t.
    107     ///
    108     /// \param k How many paths are we looking for?
    109     ///
    110     int run(int k) {
    111       int i = min_cost_flow.run(k);
    112 
    113       //Let's find the paths
    114       //We put the paths into stl vectors (as an inner representation).
    115       //In the meantime we lose the information stored in 'reversed'.
    116       //We suppose the lengths to be positive now.
    117 
    118       //We don't want to change the flow of min_cost_flow, so we make a copy
    119       //The name here suggests that the flow has only 0/1 values.
    120       EdgeIntMap reversed(G);
    121 
    122       for(typename Graph::EdgeIt e(G); e!=INVALID; ++e)
    123         reversed[e] = min_cost_flow.getFlow()[e];
    124      
     294    ///
     295    /// \param k The number of paths to be found.
     296    ///
     297    /// \return \c k if there are at least \c k edge-disjoint paths
     298    /// from \c s to \c t. Otherwise it returns the number of
     299    /// edge-disjoint paths found.
     300    ///
     301    /// \note Apart from the return value, <tt>s.run(k)</tt> is just a
     302    /// shortcut of the following code.
     303    /// \code
     304    ///   s.init();
     305    ///   s.findFlow(k);
     306    ///   s.findPaths();
     307    /// \endcode
     308    int run(int k = 2) {
     309      init();
     310      findFlow(k);
     311      findPaths();
     312      return _path_num;
     313    }
     314
     315    /// \brief Initializes the algorithm.
     316    ///
     317    /// Initializes the algorithm.
     318    void init() {
     319      // Initializing maps
     320      if (!_flow) {
     321        _flow = new FlowMap(_graph);
     322        _local_flow = true;
     323      }
     324      if (!_potential) {
     325        _potential = new PotentialMap(_graph);
     326        _local_potential = true;
     327      }
     328      for (EdgeIt e(_graph); e != INVALID; ++e) (*_flow)[e] = 0;
     329      for (NodeIt n(_graph); n != INVALID; ++n) (*_potential)[n] = 0;
     330
     331      _dijkstra = new ResidualDijkstra( _graph, *_flow, _length,
     332                                        *_potential, _pred,
     333                                        _source, _target );
     334    }
     335
     336    /// \brief Executes the successive shortest path algorithm to find
     337    /// an optimal flow.
     338    ///
     339    /// Executes the successive shortest path algorithm to find a
     340    /// minimum cost flow, which is the union of \c k or less
     341    /// edge-disjoint paths.
     342    ///
     343    /// \return \c k if there are at least \c k edge-disjoint paths
     344    /// from \c s to \c t. Otherwise it returns the number of
     345    /// edge-disjoint paths found.
     346    ///
     347    /// \pre \ref init() must be called before using this function.
     348    int findFlow(int k = 2) {
     349      // Finding shortest paths
     350      _path_num = 0;
     351      while (_path_num < k) {
     352        // Running Dijkstra
     353        if (!_dijkstra->run()) break;
     354        ++_path_num;
     355
     356        // Setting the flow along the found shortest path
     357        Node u = _target;
     358        Edge e;
     359        while ((e = _pred[u]) != INVALID) {
     360          if (u == _graph.target(e)) {
     361            (*_flow)[e] = 1;
     362            u = _graph.source(e);
     363          } else {
     364            (*_flow)[e] = 0;
     365            u = _graph.target(e);
     366          }
     367        }
     368      }
     369      return _path_num;
     370    }
     371   
     372    /// \brief Computes the paths from the flow.
     373    ///
     374    /// Computes the paths from the flow.
     375    ///
     376    /// \pre \ref init() and \ref findFlow() must be called before using
     377    /// this function.
     378    void findPaths() {
     379      // Creating the residual flow map (the union of the paths not
     380      // found so far)
     381      FlowMap res_flow(*_flow);
     382
    125383      paths.clear();
    126       paths.resize(k);
    127       for (int j=0; j<i; ++j){
    128         Node n=s;
    129 
    130         while (n!=t){
    131 
    132           OutEdgeIt e(G, n);
    133          
    134           while (!reversed[e]){
    135             ++e;
    136           }
    137           n = G.target(e);
    138           paths[j].addBack(e);
    139           reversed[e] = 1-reversed[e];
    140         }
    141        
    142       }
    143       return i;
    144     }
    145 
    146    
    147     /// \brief Returns the total length of the paths.
    148     ///
    149     /// This function gives back the total length of the found paths.
    150     Length totalLength(){
    151       return min_cost_flow.totalLength();
    152     }
    153 
    154     /// \brief Returns the found flow.
    155     ///
    156     /// This function returns a const reference to the EdgeMap \c flow.
    157     const EdgeIntMap &getFlow() const { return min_cost_flow.flow;}
    158 
    159     /// \brief Returns the optimal dual solution
    160     ///
    161     /// This function returns a const reference to the NodeMap \c
    162     /// potential (the dual solution).
    163     const EdgeIntMap &getPotential() const { return min_cost_flow.potential;}
    164 
    165     /// \brief Checks whether the complementary slackness holds.
    166     ///
    167     /// This function checks, whether the given solution is optimal.
    168     /// Currently this function only checks optimality, doesn't bother
    169     /// with feasibility.  It is meant for testing purposes.
    170     bool checkComplementarySlackness(){
    171       return min_cost_flow.checkComplementarySlackness();
    172     }
    173 
    174     typedef SimplePath<Graph> Path;
    175 
    176     /// \brief Read the found paths.
    177     ///
    178     /// This function gives back the \c j-th path in argument p.
    179     /// Assumes that \c run() has been run and nothing has changed
    180     /// since then.
    181     ///
    182     /// \warning It is assumed that \c p is constructed to be a path
    183     /// of graph \c G.  If \c j is not less than the result of
    184     /// previous \c run, then the result here will be an empty path
    185     /// (\c j can be 0 as well).
    186     ///
    187     /// \param j Which path you want to get from the found paths (in a
    188     /// real application you would get the found paths iteratively).
    189     Path path(int j) const {
    190       return paths[j];
    191     }
    192 
    193     /// \brief Gives back the number of the paths.
    194     ///
    195     /// Gives back the number of the constructed paths.
     384      paths.resize(_path_num);
     385      for (int i = 0; i < _path_num; ++i) {
     386        Node n = _source;
     387        while (n != _target) {
     388          OutEdgeIt e(_graph, n);
     389          for ( ; res_flow[e] == 0; ++e) ;
     390          n = _graph.target(e);
     391          paths[i].addBack(e);
     392          res_flow[e] = 0;
     393        }
     394      }
     395    }
     396
     397    /// @}
     398
     399    /// \name Query Functions
     400    /// The result of the algorithm can be obtained using these
     401    /// functions.
     402    /// \n The algorithm should be executed before using them.
     403
     404    /// @{
     405
     406    /// \brief Returns a const reference to the edge map storing the
     407    /// found flow.
     408    ///
     409    /// Returns a const reference to the edge map storing the flow that
     410    /// is the union of the found edge-disjoint paths.
     411    ///
     412    /// \pre \ref run() or findFlow() must be called before using this
     413    /// function.
     414    const FlowMap& flowMap() const {
     415      return *_flow;
     416    }
     417
     418    /// \brief Returns a const reference to the node map storing the
     419    /// found potentials (the dual solution).
     420    ///
     421    /// Returns a const reference to the node map storing the found
     422    /// potentials that provide the dual solution of the underlying
     423    /// minimum cost flow problem.
     424    ///
     425    /// \pre \ref run() or findFlow() must be called before using this
     426    /// function.
     427    const PotentialMap& potentialMap() const {
     428      return *_potential;
     429    }
     430
     431    /// \brief Returns the flow on the given edge.
     432    ///
     433    /// Returns the flow on the given edge.
     434    /// It is \c 1 if the edge is involved in one of the found paths,
     435    /// otherwise it is \c 0.
     436    ///
     437    /// \pre \ref run() or findFlow() must be called before using this
     438    /// function.
     439    int flow(const Edge& edge) const {
     440      return (*_flow)[edge];
     441    }
     442
     443    /// \brief Returns the potential of the given node.
     444    ///
     445    /// Returns the potential of the given node.
     446    ///
     447    /// \pre \ref run() or findFlow() must be called before using this
     448    /// function.
     449    Length potential(const Node& node) const {
     450      return (*_potential)[node];
     451    }
     452
     453    /// \brief Returns the total length (cost) of the found paths (flow).
     454    ///
     455    /// Returns the total length (cost) of the found paths (flow).
     456    /// The complexity of the function is \f$ O(e) \f$.
     457    ///
     458    /// \pre \ref run() or findFlow() must be called before using this
     459    /// function.
     460    Length totalLength() const {
     461      Length c = 0;
     462      for (EdgeIt e(_graph); e != INVALID; ++e)
     463        c += (*_flow)[e] * _length[e];
     464      return c;
     465    }
     466
     467    /// \brief Returns the number of the found paths.
     468    ///
     469    /// Returns the number of the found paths.
     470    ///
     471    /// \pre \ref run() or findFlow() must be called before using this
     472    /// function.
    196473    int pathNum() const {
    197       return paths.size();
    198     }
     474      return _path_num;
     475    }
     476
     477    /// \brief Returns a const reference to the specified path.
     478    ///
     479    /// Returns a const reference to the specified path.
     480    ///
     481    /// \param i The function returns the \c i-th path.
     482    /// \c i must be between \c 0 and <tt>%pathNum()-1</tt>.
     483    ///
     484    /// \pre \ref run() or findPaths() must be called before using this
     485    /// function.
     486    Path path(int i) const {
     487      return paths[i];
     488    }
     489
     490    /// @}
    199491
    200492  }; //class Suurballe
  • test/suurballe_test.cc

    r2553 r2586  
    1818
    1919#include <iostream>
     20#include <fstream>
     21
    2022#include <lemon/list_graph.h>
     23#include <lemon/graph_reader.h>
     24#include <lemon/path.h>
    2125#include <lemon/suurballe.h>
    22 //#include <path.h>
     26
    2327#include "test_tools.h"
    2428
    2529using namespace lemon;
    2630
     31// Checks the feasibility of the flow
     32template <typename Graph, typename FlowMap>
     33bool checkFlow( const Graph& gr, const FlowMap& flow,
     34                typename Graph::Node s, typename Graph::Node t,
     35                int value )
     36{
     37  GRAPH_TYPEDEFS(typename Graph);
     38  for (EdgeIt e(gr); e != INVALID; ++e)
     39    if (!(flow[e] == 0 || flow[e] == 1)) return false;
    2740
    28 bool passed = true;
     41  for (NodeIt n(gr); n != INVALID; ++n) {
     42    int sum = 0;
     43    for (OutEdgeIt e(gr, n); e != INVALID; ++e)
     44      sum += flow[e];
     45    for (InEdgeIt e(gr, n); e != INVALID; ++e)
     46      sum -= flow[e];
     47    if (n == s && sum != value) return false;
     48    if (n == t && sum != -value) return false;
     49    if (n != s && n != t && sum != 0) return false;
     50  }
     51
     52  return true;
     53}
     54
     55// Checks the optimalitiy of the flow
     56template < typename Graph, typename CostMap,
     57           typename FlowMap, typename PotentialMap >
     58bool checkOptimality( const Graph& gr, const CostMap& cost,
     59                      const FlowMap& flow, const PotentialMap& pi )
     60{
     61  // Checking the Complementary Slackness optimality condition
     62  GRAPH_TYPEDEFS(typename Graph);
     63  bool opt = true;
     64  for (EdgeIt e(gr); e != INVALID; ++e) {
     65    typename CostMap::Value red_cost =
     66      cost[e] + pi[gr.source(e)] - pi[gr.target(e)];
     67    opt = (flow[e] == 0 && red_cost >= 0) ||
     68          (flow[e] == 1 && red_cost <= 0);
     69    if (!opt) break;
     70  }
     71  return opt;
     72}
     73
     74// Checks a path
     75template < typename Graph, typename Path >
     76bool checkPath( const Graph& gr, const Path& path,
     77                typename Graph::Node s, typename Graph::Node t)
     78{
     79  // Checking the Complementary Slackness optimality condition
     80  GRAPH_TYPEDEFS(typename Graph);
     81  Node n = s;
     82  for (int i = 0; i < path.length(); ++i) {
     83    if (gr.source(path.nth(i)) != n) return false;
     84    n = gr.target(path.nth(i));
     85  }
     86  return n == t;
     87}
    2988
    3089
    3190int main()
    3291{
    33   typedef ListGraph Graph;
    34   typedef Graph::Node Node;
    35   typedef Graph::Edge Edge;
     92  GRAPH_TYPEDEFS(ListGraph);
    3693
    37   Graph graph;
     94  // Reading the test graph
     95  ListGraph graph;
     96  ListGraph::EdgeMap<int> length(graph);
     97  Node source, target;
    3898
    39   //Ahuja könyv példája
     99  std::string fname;
     100  if(getenv("srcdir"))
     101    fname = std::string(getenv("srcdir"));
     102  else fname = ".";
     103  fname += "/test/min_cost_flow_test.lgf";
    40104
    41   Node s=graph.addNode();
    42   Node v1=graph.addNode(); 
    43   Node v2=graph.addNode();
    44   Node v3=graph.addNode();
    45   Node v4=graph.addNode();
    46   Node v5=graph.addNode();
    47   Node t=graph.addNode();
     105  std::ifstream input(fname.c_str());
     106  check(input, "Input file '" << fname << "' not found");
     107  GraphReader<ListGraph>(input, graph).
     108    readEdgeMap("cost", length).
     109    readNode("source", source).
     110    readNode("target", target).
     111    run();
     112  input.close();
     113 
     114  // Finding 2 paths
     115  {
     116    Suurballe<ListGraph> suurballe(graph, length, source, target);
     117    check(suurballe.run(2) == 2, "Wrong number of paths");
     118    check(checkFlow(graph, suurballe.flowMap(), source, target, 2),
     119          "The flow is not feasible");
     120    check(suurballe.totalLength() == 510, "The flow is not optimal");
     121    check(checkOptimality(graph, length, suurballe.flowMap(),
     122                          suurballe.potentialMap()),
     123          "Wrong potentials");
     124    for (int i = 0; i < suurballe.pathNum(); ++i)
     125      check(checkPath(graph, suurballe.path(i), source, target),
     126            "Wrong path");
     127  }
    48128
    49   Edge s_v1=graph.addEdge(s, v1);
    50   Edge v1_v2=graph.addEdge(v1, v2);
    51   Edge s_v3=graph.addEdge(s, v3);
    52   Edge v2_v4=graph.addEdge(v2, v4);
    53   Edge v2_v5=graph.addEdge(v2, v5);
    54   Edge v3_v5=graph.addEdge(v3, v5);
    55   Edge v4_t=graph.addEdge(v4, t);
    56   Edge v5_t=graph.addEdge(v5, t);
    57  
     129  // Finding 3 paths
     130  {
     131    Suurballe<ListGraph> suurballe(graph, length, source, target);
     132    check(suurballe.run(3) == 3, "Wrong number of paths");
     133    check(checkFlow(graph, suurballe.flowMap(), source, target, 3),
     134          "The flow is not feasible");
     135    check(suurballe.totalLength() == 1040, "The flow is not optimal");
     136    check(checkOptimality(graph, length, suurballe.flowMap(),
     137                          suurballe.potentialMap()),
     138          "Wrong potentials");
     139    for (int i = 0; i < suurballe.pathNum(); ++i)
     140      check(checkPath(graph, suurballe.path(i), source, target),
     141            "Wrong path");
     142  }
    58143
    59   Graph::EdgeMap<int> length(graph);
     144  // Finding 5 paths (only 3 can be found)
     145  {
     146    Suurballe<ListGraph> suurballe(graph, length, source, target);
     147    check(suurballe.run(5) == 3, "Wrong number of paths");
     148    check(checkFlow(graph, suurballe.flowMap(), source, target, 3),
     149          "The flow is not feasible");
     150    check(suurballe.totalLength() == 1040, "The flow is not optimal");
     151    check(checkOptimality(graph, length, suurballe.flowMap(),
     152                          suurballe.potentialMap()),
     153          "Wrong potentials");
     154    for (int i = 0; i < suurballe.pathNum(); ++i)
     155      check(checkPath(graph, suurballe.path(i), source, target),
     156            "Wrong path");
     157  }
    60158
    61   length.set(s_v1, 6);
    62   length.set(v1_v2, 4);
    63   length.set(s_v3, 10);
    64   length.set(v2_v4, 5);
    65   length.set(v2_v5, 1);
    66   length.set(v3_v5, 5);
    67   length.set(v4_t, 8);
    68   length.set(v5_t, 8);
    69 
    70   std::cout << "Minlengthpaths algorithm test..." << std::endl;
    71 
    72  
    73   int k=3;
    74   Suurballe< Graph, Graph::EdgeMap<int> >
    75     surb_test(graph, length, s, t);
    76 
    77   check(  surb_test.run(k) == 2 && surb_test.totalLength() == 46,
    78           "Two paths, total length should be 46");
    79 
    80   check(  surb_test.checkComplementarySlackness(),
    81           "Complementary slackness conditions are not met.");
    82 
    83   //  typedef DirPath<Graph> DPath;
    84   //  DPath P(graph);
    85 
    86   /*
    87   surb_test.getPath(P,0);
    88   check(P.length() == 4, "First path should contain 4 edges."); 
    89   std::cout<<P.length()<<std::endl;
    90   surb_test.getPath(P,1);
    91   check(P.length() == 3, "Second path: 3 edges.");
    92   std::cout<<P.length()<<std::endl;
    93   */ 
    94 
    95   k=1;
    96   check(  surb_test.run(k) == 1 && surb_test.totalLength() == 19,
    97           "One path, total length should be 19");
    98 
    99   check(  surb_test.checkComplementarySlackness(),
    100           "Complementary slackness conditions are not met.");
    101  
    102   //  surb_test.getPath(P,0);
    103   //  check(P.length() == 4, "First path should contain 4 edges."); 
    104 
    105   std::cout << (passed ? "All tests passed." : "Some of the tests failed!!!")
    106             << std::endl;
    107 
    108   return passed ? 0 : 1;
    109 
     159  return 0;
    110160}
Note: See TracChangeset for help on using the changeset viewer.