lemon/min_cost_arborescence.h
changeset 2020 332245399dc6
child 2025 93fcadf94ab0
equal deleted inserted replaced
-1:000000000000 0:88bf8b76d70e
       
     1 /* -*- C++ -*-
       
     2  *
       
     3  * This file is a part of LEMON, a generic C++ optimization library
       
     4  *
       
     5  * Copyright (C) 2003-2006
       
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
       
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
       
     8  *
       
     9  * Permission to use, modify and distribute this software is granted
       
    10  * provided that this copyright notice appears in all copies. For
       
    11  * precise terms see the accompanying LICENSE file.
       
    12  *
       
    13  * This software is provided "AS IS" with no warranty of any kind,
       
    14  * express or implied, and with no claim as to its suitability for any
       
    15  * purpose.
       
    16  *
       
    17  */
       
    18 
       
    19 #ifndef LEMON_MIN_COST_ARBORESCENCE_H
       
    20 #define LEMON_MIN_COST_ARBORESCENCE_H
       
    21 
       
    22 ///\ingroup spantree
       
    23 ///\file
       
    24 ///\brief Minimum Cost Arborescence algorithm.
       
    25 
       
    26 #include <vector>
       
    27 
       
    28 #include <lemon/list_graph.h>
       
    29 
       
    30 namespace lemon {
       
    31 
       
    32 
       
    33   /// \brief Default traits class of MinCostArborescence class.
       
    34   /// 
       
    35   /// Default traits class of MinCostArborescence class.
       
    36   /// \param _Graph Graph type.
       
    37   /// \param _CostMap Type of cost map.
       
    38   template <class _Graph, class _CostMap>
       
    39   struct MinCostArborescenceDefaultTraits{
       
    40 
       
    41     /// \brief The graph type the algorithm runs on. 
       
    42     typedef _Graph Graph;
       
    43 
       
    44     /// \brief The type of the map that stores the edge costs.
       
    45     ///
       
    46     /// The type of the map that stores the edge costs.
       
    47     /// It must meet the \ref concept::ReadMap "ReadMap" concept.
       
    48     typedef _CostMap CostMap;
       
    49 
       
    50     /// \brief The value type of the costs.
       
    51     ///
       
    52     /// The value type of the costs.
       
    53     typedef typename CostMap::Value Value;
       
    54 
       
    55     /// \brief The type of the map that stores which edges are 
       
    56     /// in the arborescence.
       
    57     ///
       
    58     /// The type of the map that stores which edges are in the arborescence.
       
    59     /// It must meet the \ref concept::ReadWriteMap "ReadWriteMap" concept.
       
    60     /// Initially it will be setted to false on each edge. The algorithm
       
    61     /// may set each value one time to true and maybe after it to false again.
       
    62     /// Therefore you cannot use maps like BackInserteBoolMap with this
       
    63     /// algorithm.   
       
    64     typedef typename Graph::template EdgeMap<bool> ArborescenceMap; 
       
    65 
       
    66     /// \brief Instantiates a ArborescenceMap.
       
    67     ///
       
    68     /// This function instantiates a \ref ArborescenceMap. 
       
    69     /// \param _graph is the graph, to which we would like to define the 
       
    70     /// ArborescenceMap.
       
    71     static ArborescenceMap *createArborescenceMap(const Graph &_graph){
       
    72       return new ArborescenceMap(_graph);
       
    73     }
       
    74 
       
    75   };
       
    76 
       
    77   /// \ingroup spantree
       
    78   ///
       
    79   /// \brief %MinCostArborescence algorithm class.
       
    80   ///
       
    81   /// This class provides an efficient implementation of 
       
    82   /// %MinCostArborescence algorithm. The arborescence is a tree 
       
    83   /// which is directed from a given source node of the graph. One or
       
    84   /// more sources should be given for the algorithm and it will calculate
       
    85   /// the minimum cost subgraph which are union of arborescences with the
       
    86   /// given sources and spans all the nodes which are reachable from the
       
    87   /// sources. The time complexity of the algorithm is O(n^2 + e).
       
    88   ///
       
    89   /// \param _Graph The graph type the algorithm runs on. The default value
       
    90   /// is \ref ListGraph. The value of _Graph is not used directly by
       
    91   /// MinCostArborescence, it is only passed to 
       
    92   /// \ref MinCostArborescenceDefaultTraits.
       
    93   /// \param _CostMap This read-only EdgeMap determines the costs of the
       
    94   /// edges. It is read once for each edge, so the map may involve in
       
    95   /// relatively time consuming process to compute the edge cost if
       
    96   /// it is necessary. The default map type is \ref
       
    97   /// concept::StaticGraph::EdgeMap "Graph::EdgeMap<int>".  The value
       
    98   /// of _CostMap is not used directly by MinCostArborescence, 
       
    99   /// it is only passed to \ref MinCostArborescenceDefaultTraits.  
       
   100   /// \param _Traits Traits class to set various data types used 
       
   101   /// by the algorithm.  The default traits class is 
       
   102   /// \ref MinCostArborescenceDefaultTraits
       
   103   /// "MinCostArborescenceDefaultTraits<_Graph,_CostMap>".  See \ref
       
   104   /// MinCostArborescenceDefaultTraits for the documentation of a 
       
   105   /// MinCostArborescence traits class.
       
   106   ///
       
   107   /// \author Balazs Dezso
       
   108 #ifndef DOXYGEN
       
   109   template <typename _Graph = ListGraph, 
       
   110             typename _CostMap = typename _Graph::template EdgeMap<int>,
       
   111             typename _Traits = 
       
   112             MinCostArborescenceDefaultTraits<_Graph, _CostMap> >
       
   113 #else 
       
   114   template <typename _Graph, typename _CostMap, typedef _Traits>
       
   115 #endif
       
   116   class MinCostArborescence {
       
   117   public:
       
   118 
       
   119     /// \brief \ref Exception for uninitialized parameters.
       
   120     ///
       
   121     /// This error represents problems in the initialization
       
   122     /// of the parameters of the algorithms.    
       
   123     class UninitializedParameter : public lemon::UninitializedParameter {
       
   124     public:
       
   125       virtual const char* exceptionName() const {
       
   126 	return "lemon::MinCostArborescence::UninitializedParameter";
       
   127       }
       
   128     };
       
   129 
       
   130     /// The traits.
       
   131     typedef _Traits Traits;
       
   132     /// The type of the underlying graph.
       
   133     typedef typename Traits::Graph Graph;
       
   134     /// The type of the map that stores the edge costs.
       
   135     typedef typename Traits::CostMap CostMap;
       
   136     ///The type of the costs of the edges.
       
   137     typedef typename Traits::Value Value;
       
   138     ///The type of the map that stores which edges are in the arborescence.
       
   139     typedef typename Traits::ArborescenceMap ArborescenceMap;
       
   140 
       
   141   protected:
       
   142 
       
   143     typedef typename Graph::Node Node;
       
   144     typedef typename Graph::Edge Edge;
       
   145     typedef typename Graph::NodeIt NodeIt;
       
   146     typedef typename Graph::EdgeIt EdgeIt;
       
   147     typedef typename Graph::InEdgeIt InEdgeIt;
       
   148     typedef typename Graph::OutEdgeIt OutEdgeIt;
       
   149 
       
   150     struct CostEdge {
       
   151 
       
   152       Edge edge;
       
   153       Value value;
       
   154 
       
   155       CostEdge() {}
       
   156       CostEdge(Edge _edge, Value _value) : edge(_edge), value(_value) {}
       
   157 
       
   158     };
       
   159 
       
   160     const Graph* graph;
       
   161     const CostMap* cost;
       
   162 
       
   163     ArborescenceMap* _arborescence_map;
       
   164     bool local_arborescence_map;
       
   165 
       
   166     typedef typename Graph::template NodeMap<int> LevelMap;
       
   167     LevelMap *_level;
       
   168 
       
   169     typedef typename Graph::template NodeMap<CostEdge> CostEdgeMap;
       
   170     CostEdgeMap *_cost_edges; 
       
   171 
       
   172     struct StackLevel {
       
   173 
       
   174       std::vector<CostEdge> edges;
       
   175       int node_level;
       
   176 
       
   177     };
       
   178 
       
   179     std::vector<StackLevel> level_stack;    
       
   180     std::vector<Node> queue;
       
   181 
       
   182     int node_counter;
       
   183 
       
   184   public:
       
   185 
       
   186     /// \name Named template parameters
       
   187 
       
   188     /// @{
       
   189 
       
   190     template <class T>
       
   191     struct DefArborescenceMapTraits : public Traits {
       
   192       typedef T ArborescenceMap;
       
   193       static ArborescenceMap *createArborescenceMap(const Graph &)
       
   194       {
       
   195 	throw UninitializedParameter();
       
   196       }
       
   197     };
       
   198 
       
   199     /// \brief \ref named-templ-param "Named parameter" for 
       
   200     /// setting ArborescenceMap type
       
   201     ///
       
   202     /// \ref named-templ-param "Named parameter" for setting 
       
   203     /// ArborescenceMap type
       
   204     template <class T>
       
   205     struct DefArborescenceMap 
       
   206       : public MinCostArborescence<Graph, CostMap,
       
   207                                    DefArborescenceMapTraits<T> > {
       
   208       typedef MinCostArborescence<Graph, CostMap, 
       
   209                                    DefArborescenceMapTraits<T> > Create;
       
   210     };
       
   211     
       
   212     /// @}
       
   213 
       
   214     /// \brief Constructor.
       
   215     ///
       
   216     /// \param _graph The graph the algorithm will run on.
       
   217     /// \param _cost The cost map used by the algorithm.
       
   218     MinCostArborescence(const Graph& _graph, const CostMap& _cost) 
       
   219       : graph(&_graph), cost(&_cost),
       
   220         _arborescence_map(0), local_arborescence_map(false), 
       
   221         _level(0), _cost_edges(0) {}
       
   222 
       
   223     /// \brief Destructor.
       
   224     ~MinCostArborescence() {
       
   225       destroyStructures();
       
   226     }
       
   227 
       
   228     /// \brief Sets the arborescence map.
       
   229     /// 
       
   230     /// Sets the arborescence map.
       
   231     /// \return \c (*this)
       
   232     MinCostArborescence& arborescenceMap(ArborescenceMap& m) {
       
   233       _arborescence_map = &m;
       
   234       return *this;
       
   235     }
       
   236 
       
   237     /// \name Query Functions
       
   238     /// The result of the %MinCostArborescence algorithm can be obtained 
       
   239     /// using these functions.\n
       
   240     /// Before the use of these functions,
       
   241     /// either run() or start() must be called.
       
   242     
       
   243     /// @{
       
   244 
       
   245     /// \brief Returns a reference to the arborescence map.
       
   246     ///
       
   247     /// Returns a reference to the arborescence map.
       
   248     const ArborescenceMap& arborescenceMap() const {
       
   249       return *_arborescence_map;
       
   250     }
       
   251 
       
   252     /// \brief Returns true if the edge is in the arborescence.
       
   253     ///
       
   254     /// Returns true if the edge is in the arborescence.
       
   255     /// \param edge The edge of the graph.
       
   256     /// \pre \ref run() must be called before using this function.
       
   257     bool arborescenceEdge(Edge edge) const {
       
   258       return (*_arborescence_map)[edge];
       
   259     }
       
   260  
       
   261     /// \brief Returns the cost of the arborescence.
       
   262     ///
       
   263     /// Returns the cost of the arborescence.
       
   264    Value arborescenceCost() const {
       
   265       Value sum = 0;
       
   266       for (EdgeIt it(*graph); it != INVALID; ++it) {
       
   267         if (arborescenceEdge(it)) {
       
   268           sum += (*cost)[it];
       
   269         }
       
   270       }
       
   271       return sum;
       
   272     }
       
   273 
       
   274     /// @}
       
   275     
       
   276     /// \name Execution control
       
   277     /// The simplest way to execute the algorithm is to use
       
   278     /// one of the member functions called \c run(...). \n
       
   279     /// If you need more control on the execution,
       
   280     /// first you must call \ref init(), then you can add several 
       
   281     /// source nodes with \ref addSource().
       
   282     /// Finally \ref start() will perform the actual path
       
   283     /// computation.
       
   284 
       
   285     ///@{
       
   286 
       
   287     /// \brief Initializes the internal data structures.
       
   288     ///
       
   289     /// Initializes the internal data structures.
       
   290     ///
       
   291     void init() {
       
   292       initStructures();
       
   293       for (NodeIt it(*graph); it != INVALID; ++it) {
       
   294         (*_cost_edges)[it].edge = INVALID;
       
   295         (*_level)[it] = -3; 
       
   296       }
       
   297       for (EdgeIt it(*graph); it != INVALID; ++it) {
       
   298         _arborescence_map->set(it, false);
       
   299       }
       
   300     }
       
   301 
       
   302     /// \brief Adds a new source node.
       
   303     ///
       
   304     /// Adds a new source node to the algorithm.
       
   305     void addSource(Node source) {
       
   306       std::vector<Node> nodes;
       
   307       nodes.push_back(source);
       
   308       while (!nodes.empty()) {
       
   309         Node node = nodes.back();
       
   310         nodes.pop_back();
       
   311         for (OutEdgeIt it(*graph, node); it != INVALID; ++it) {
       
   312           if ((*_level)[graph->target(it)] == -3) {
       
   313             (*_level)[graph->target(it)] = -2;
       
   314             nodes.push_back(graph->target(it));
       
   315             queue.push_back(graph->target(it));
       
   316           }
       
   317         }
       
   318       }
       
   319       (*_level)[source] = -1;
       
   320     }
       
   321 
       
   322     /// \brief Processes the next node in the priority queue.
       
   323     ///
       
   324     /// Processes the next node in the priority queue.
       
   325     ///
       
   326     /// \return The processed node.
       
   327     ///
       
   328     /// \warning The queue must not be empty!
       
   329     Node processNextNode() {
       
   330       node_counter = 0;
       
   331       Node node = queue.back();
       
   332       queue.pop_back();
       
   333       if ((*_level)[node] == -2) {
       
   334         Edge edge = prepare(node);
       
   335         while ((*_level)[graph->source(edge)] != -1) {
       
   336           if ((*_level)[graph->source(edge)] >= 0) {
       
   337             edge = contract(bottom((*_level)[graph->source(edge)]));
       
   338           } else {
       
   339             edge = prepare(graph->source(edge));
       
   340           }
       
   341         }
       
   342         finalize(graph->target(edge));
       
   343         level_stack.clear();        
       
   344       }
       
   345       return node;
       
   346     }
       
   347 
       
   348     /// \brief Returns the number of the nodes to be processed.
       
   349     ///
       
   350     /// Returns the number of the nodes to be processed.
       
   351     int queueSize() const {
       
   352       return queue.size();
       
   353     }
       
   354 
       
   355     /// \brief Returns \c false if there are nodes to be processed.
       
   356     ///
       
   357     /// Returns \c false if there are nodes to be processed.
       
   358     bool emptyQueue() const {
       
   359       return queue.empty();
       
   360     }
       
   361 
       
   362     /// \brief Executes the algorithm.
       
   363     ///
       
   364     /// Executes the algorithm.
       
   365     ///
       
   366     /// \pre init() must be called and at least one node should be added
       
   367     /// with addSource() before using this function.
       
   368     ///
       
   369     ///\note mca.start() is just a shortcut of the following code.
       
   370     ///\code
       
   371     ///while (!mca.emptyQueue()) {
       
   372     ///  mca.processNextNode();
       
   373     ///}
       
   374     ///\endcode
       
   375     void start() {
       
   376       while (!emptyQueue()) {
       
   377         processNextNode();
       
   378       }
       
   379     }
       
   380 
       
   381     /// \brief Runs %MinCostArborescence algorithm from node \c s.
       
   382     /// 
       
   383     /// This method runs the %MinCostArborescence algorithm from 
       
   384     /// a root node \c s.
       
   385     ///
       
   386     ///\note mca.run(s) is just a shortcut of the following code.
       
   387     ///\code
       
   388     ///mca.init();
       
   389     ///mca.addSource(s);
       
   390     ///mca.start();
       
   391     ///\endcode
       
   392     void run(Node node) {
       
   393       init();
       
   394       addSource(node);
       
   395       start();
       
   396     }
       
   397 
       
   398     ///@}
       
   399 
       
   400   protected:
       
   401 
       
   402     void initStructures() {
       
   403       if (!_arborescence_map) {
       
   404         local_arborescence_map = true;
       
   405         _arborescence_map = Traits::createArborescenceMap(*graph);
       
   406       }
       
   407       if (!_level) {
       
   408         _level = new LevelMap(*graph);
       
   409       }
       
   410       if (!_cost_edges) {
       
   411         _cost_edges = new CostEdgeMap(*graph);
       
   412       }
       
   413     }
       
   414 
       
   415     void destroyStructures() {
       
   416       if (_level) {
       
   417         delete _level;
       
   418       }
       
   419       if (!_cost_edges) {
       
   420         delete _cost_edges;
       
   421       }
       
   422       if (local_arborescence_map) {
       
   423         delete _arborescence_map;
       
   424       }
       
   425     }
       
   426 
       
   427     Edge prepare(Node node) {
       
   428       std::vector<Node> nodes;
       
   429       (*_level)[node] = node_counter;
       
   430       for (InEdgeIt it(*graph, node); it != INVALID; ++it) {
       
   431         Edge edge = it;
       
   432         Value value = (*cost)[it];
       
   433         if (graph->source(edge) == node || 
       
   434             (*_level)[graph->source(edge)] == -3) continue;
       
   435         if ((*_cost_edges)[graph->source(edge)].edge == INVALID) {
       
   436           (*_cost_edges)[graph->source(edge)].edge = edge;
       
   437           (*_cost_edges)[graph->source(edge)].value = value;
       
   438           nodes.push_back(graph->source(edge));
       
   439         } else {
       
   440           if ((*_cost_edges)[graph->source(edge)].value > value) {
       
   441             (*_cost_edges)[graph->source(edge)].edge = edge;
       
   442             (*_cost_edges)[graph->source(edge)].value = value;
       
   443           }
       
   444         }      
       
   445       }
       
   446       CostEdge minimum = (*_cost_edges)[nodes[0]]; 
       
   447       for (int i = 1; i < (int)nodes.size(); ++i) {
       
   448         if ((*_cost_edges)[nodes[i]].value < minimum.value) {
       
   449           minimum = (*_cost_edges)[nodes[i]];
       
   450         }
       
   451       }
       
   452       StackLevel level;
       
   453       level.node_level = node_counter;
       
   454       for (int i = 0; i < (int)nodes.size(); ++i) {
       
   455         (*_cost_edges)[nodes[i]].value -= minimum.value;
       
   456         level.edges.push_back((*_cost_edges)[nodes[i]]);
       
   457         (*_cost_edges)[nodes[i]].edge = INVALID;
       
   458       }
       
   459       level_stack.push_back(level);
       
   460       ++node_counter;
       
   461       _arborescence_map->set(minimum.edge, true);
       
   462       return minimum.edge;
       
   463     }
       
   464   
       
   465     Edge contract(int node_bottom) {
       
   466       std::vector<Node> nodes;
       
   467       while (!level_stack.empty() && 
       
   468              level_stack.back().node_level >= node_bottom) {
       
   469         for (int i = 0; i < (int)level_stack.back().edges.size(); ++i) {
       
   470           Edge edge = level_stack.back().edges[i].edge;
       
   471           Value value = level_stack.back().edges[i].value;
       
   472           if ((*_level)[graph->source(edge)] >= node_bottom) continue;
       
   473           if ((*_cost_edges)[graph->source(edge)].edge == INVALID) {
       
   474             (*_cost_edges)[graph->source(edge)].edge = edge;
       
   475             (*_cost_edges)[graph->source(edge)].value = value;
       
   476             nodes.push_back(graph->source(edge));
       
   477           } else {
       
   478             if ((*_cost_edges)[graph->source(edge)].value > value) {
       
   479               (*_cost_edges)[graph->source(edge)].edge = edge;
       
   480               (*_cost_edges)[graph->source(edge)].value = value;
       
   481             }
       
   482           }
       
   483         }
       
   484         level_stack.pop_back();
       
   485       }
       
   486       CostEdge minimum = (*_cost_edges)[nodes[0]]; 
       
   487       for (int i = 1; i < (int)nodes.size(); ++i) {
       
   488         if ((*_cost_edges)[nodes[i]].value < minimum.value) {
       
   489           minimum = (*_cost_edges)[nodes[i]];
       
   490         }
       
   491       }
       
   492       StackLevel level;
       
   493       level.node_level = node_bottom;
       
   494       for (int i = 0; i < (int)nodes.size(); ++i) {
       
   495         (*_cost_edges)[nodes[i]].value -= minimum.value;
       
   496         level.edges.push_back((*_cost_edges)[nodes[i]]);
       
   497         (*_cost_edges)[nodes[i]].edge = INVALID;
       
   498       }
       
   499       level_stack.push_back(level);
       
   500       _arborescence_map->set(minimum.edge, true);
       
   501       return minimum.edge;
       
   502     }
       
   503 
       
   504     int bottom(int level) {
       
   505       int k = level_stack.size() - 1;
       
   506       while (level_stack[k].node_level > level) {
       
   507         --k;
       
   508       }
       
   509       return level_stack[k].node_level;
       
   510     }
       
   511 
       
   512     void finalize(Node source) {
       
   513       std::vector<Node> nodes;
       
   514       nodes.push_back(source);
       
   515       while (!nodes.empty()) {
       
   516         Node node = nodes.back();
       
   517         nodes.pop_back();
       
   518         for (OutEdgeIt it(*graph, node); it != INVALID; ++it) {
       
   519           if ((*_level)[graph->target(it)] >= 0 && (*_arborescence_map)[it]) {
       
   520             (*_level)[graph->target(it)] = -1;
       
   521             nodes.push_back(graph->target(it));
       
   522           } else {
       
   523             _arborescence_map->set(it, false);
       
   524           }
       
   525         }
       
   526       }
       
   527       (*_level)[source] = -1;      
       
   528     }
       
   529 
       
   530   };
       
   531 
       
   532   /// \ingroup spantree
       
   533   ///
       
   534   /// \brief Function type interface for MinCostArborescence algorithm.
       
   535   ///
       
   536   /// Function type interface for MinCostArborescence algorithm.
       
   537   /// \param graph The Graph that the algorithm runs on.
       
   538   /// \param cost The CostMap of the edges.
       
   539   /// \param source The source of the arborescence.
       
   540   /// \retval arborescence The bool EdgeMap which stores the arborescence.
       
   541   /// \return The cost of the arborescence. 
       
   542   ///
       
   543   /// \sa MinCostArborescence
       
   544   template <typename Graph, typename CostMap, typename ArborescenceMap>
       
   545   typename CostMap::Value minCostArborescence(const Graph& graph, 
       
   546                                               const CostMap& cost,
       
   547                                               typename Graph::Node source,
       
   548                                               ArborescenceMap& arborescence) {
       
   549     typename MinCostArborescence<Graph, CostMap>
       
   550       ::template DefArborescenceMap<ArborescenceMap>
       
   551       ::Create mca(graph, cost);
       
   552     mca.arborescenceMap(arborescence);
       
   553     mca.run(source);
       
   554     return mca.arborescenceCost();
       
   555   }
       
   556 
       
   557 }
       
   558 
       
   559 #endif
       
   560 
       
   561 // Hilbert - Huang