deba@2382: /* -*- C++ -*-
deba@2382:  *
deba@2382:  * This file is a part of LEMON, a generic C++ optimization library
deba@2382:  *
alpar@2391:  * Copyright (C) 2003-2007
deba@2382:  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@2382:  * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@2382:  *
deba@2382:  * Permission to use, modify and distribute this software is granted
deba@2382:  * provided that this copyright notice appears in all copies. For
deba@2382:  * precise terms see the accompanying LICENSE file.
deba@2382:  *
deba@2382:  * This software is provided "AS IS" with no warranty of any kind,
deba@2382:  * express or implied, and with no claim as to its suitability for any
deba@2382:  * purpose.
deba@2382:  *
deba@2382:  */
deba@2382: 
deba@2382: #ifndef LEMON_STEINER_H
deba@2382: #define LEMON_STEINER_H
deba@2382: 
deba@2382: ///\ingroup approx
deba@2382: ///\file
deba@2382: ///\brief Algorithm for the 2-approximation of Steiner Tree problem.
deba@2382: ///
deba@2382: 
deba@2382: #include <lemon/smart_graph.h>
deba@2382: #include <lemon/graph_utils.h>
deba@2382: #include <lemon/error.h>
deba@2382: 
deba@2382: #include <lemon/ugraph_adaptor.h>
deba@2382: #include <lemon/maps.h>
deba@2382: 
deba@2382: #include <lemon/dijkstra.h>
deba@2382: #include <lemon/prim.h>
deba@2382: 
deba@2382: 
deba@2382: namespace lemon {
deba@2382: 
deba@2382:   /// \ingroup approx
deba@2382:   
deba@2382:   /// \brief Algorithm for the 2-approximation of Steiner Tree problem
deba@2382:   ///
deba@2382:   /// The Steiner-tree problem is the next: Given a connected
deba@2382:   /// undirected graph, a cost function on the edges and a subset of
deba@2382:   /// the nodes. Construct a tree with minimum cost which covers the
deba@2382:   /// given subset of the nodes. The problem is NP-hard moreover
deba@2382:   /// it is APX-complete too.
deba@2382:   ///
deba@2382:   /// Mehlhorn's approximation algorithm is implemented in this class,
deba@2382:   /// which gives a 2-approximation for the Steiner-tree problem. The
deba@2382:   /// algorithm's time complexity is O(nlog(n)+e).
deba@2382:   template <typename UGraph,
deba@2382:             typename CostMap = typename UGraph:: template UEdgeMap<double> >
deba@2382:   class SteinerTree {
deba@2382:   public:
deba@2382:     
deba@2382:     UGRAPH_TYPEDEFS(typename UGraph)
deba@2382: 
deba@2382:     typedef typename CostMap::Value Value;
deba@2382:     
deba@2382:   private:
deba@2382: 
deba@2382:     class CompMap {
deba@2382:     public:
deba@2382:       typedef Node Key;
deba@2382:       typedef Edge Value;
deba@2382: 
deba@2382:     private:
deba@2382:       const UGraph& _graph;
deba@2382:       typename UGraph::template NodeMap<int> _comp;
deba@2382: 
deba@2382:     public:
deba@2382:       CompMap(const UGraph& graph) : _graph(graph), _comp(graph) {}
deba@2382: 
deba@2382:       void set(const Node& node, const Edge& edge) {
deba@2382:         if (edge != INVALID) {
deba@2382:           _comp.set(node, _comp[_graph.source(edge)]);
deba@2382:         } else {
deba@2382:           _comp.set(node, -1);
deba@2382:         }
deba@2382:       }
deba@2382: 
deba@2382:       int comp(const Node& node) const { return _comp[node]; }
deba@2382:       void comp(const Node& node, int value) { _comp.set(node, value); }
deba@2382:     };
deba@2382: 
deba@2382:     typedef typename UGraph::template NodeMap<Edge> PredMap;
deba@2382: 
deba@2382:     typedef ForkWriteMap<PredMap, CompMap> ForkedMap;
deba@2382: 
deba@2382: 
deba@2382:     struct External {
deba@2382:       int source, target;
deba@2382:       UEdge uedge;
deba@2382:       Value value;
deba@2382: 
deba@2382:       External(int s, int t, const UEdge& e, const Value& v)
deba@2382:         : source(s), target(t), uedge(e), value(v) {}
deba@2382:     };
deba@2382: 
deba@2382:     struct ExternalLess {
deba@2382:       bool operator()(const External& left, const External& right) const {
deba@2382:         return (left.source < right.source) || 
deba@2382:           (left.source == right.source && left.target < right.target);
deba@2382:       }
deba@2382:     };
deba@2382: 
deba@2382: 
deba@2382:     typedef typename UGraph::template NodeMap<bool> FilterMap;
deba@2382: 
deba@2382:     typedef typename UGraph::template UEdgeMap<bool> TreeMap;
deba@2382: 
deba@2382:     const UGraph& _graph;
deba@2382:     const CostMap& _cost;
deba@2382: 
deba@2382:     typename Dijkstra<UGraph, CostMap>::
deba@2382:     template DefPredMap<ForkedMap>::Create _dijkstra;
deba@2382: 
deba@2382:     PredMap* _pred;
deba@2382:     CompMap* _comp;
deba@2382:     ForkedMap* _forked;
deba@2382: 
deba@2382:     int _terminal_num;
deba@2382: 
deba@2382:     FilterMap *_filter;
deba@2382:     TreeMap *_tree;
deba@2382: 
deba@2399:     Value _value;
deba@2399: 
deba@2382:   public:
deba@2382: 
deba@2382:     /// \brief Constructor
deba@2382:     
deba@2382:     /// Constructor
deba@2382:     ///
deba@2382:     SteinerTree(const UGraph &graph, const CostMap &cost)
deba@2382:       : _graph(graph), _cost(cost), _dijkstra(graph, _cost), 
deba@2382:         _pred(0), _comp(0), _forked(0), _filter(0), _tree(0) {}
deba@2382: 
deba@2382:     /// \brief Initializes the internal data structures.
deba@2382:     ///
deba@2382:     /// Initializes the internal data structures.
deba@2382:     void init() {
deba@2382:       if (!_pred) _pred = new PredMap(_graph);
deba@2382:       if (!_comp) _comp = new CompMap(_graph);
deba@2382:       if (!_forked) _forked = new ForkedMap(*_pred, *_comp);
deba@2382:       if (!_filter) _filter = new FilterMap(_graph);
deba@2382:       if (!_tree) _tree = new TreeMap(_graph);
deba@2382:       _dijkstra.predMap(*_forked);
deba@2382:       _dijkstra.init();
deba@2382:       _terminal_num = 0;
deba@2382:       for (NodeIt it(_graph); it != INVALID; ++it) {
deba@2382:         _filter->set(it, false);
deba@2382:       }
deba@2382:     }
deba@2382: 
deba@2382:     /// \brief Adds a new terminal node.
deba@2382:     ///
deba@2382:     /// Adds a new terminal node to the Steiner-tree problem.
deba@2382:     void addTerminal(const Node& node) {
deba@2382:       if (!_dijkstra.reached(node)) {
deba@2382:         _dijkstra.addSource(node);
deba@2382:         _comp->comp(node, _terminal_num);
deba@2382:         ++_terminal_num;
deba@2382:       }
deba@2382:     }
deba@2382:     
deba@2382:     /// \brief Executes the algorithm.
deba@2382:     ///
deba@2382:     /// Executes the algorithm.
deba@2382:     ///
deba@2382:     /// \pre init() must be called and at least some nodes should be
deba@2382:     /// added with addTerminal() before using this function.
deba@2382:     ///
deba@2382:     /// This method constructs an approximation of the Steiner-Tree.
deba@2382:     void start() {
deba@2382:       _dijkstra.start();
deba@2382:       
deba@2382:       std::vector<External> externals;
deba@2382:       for (UEdgeIt it(_graph); it != INVALID; ++it) {
deba@2382:         Node s = _graph.source(it);
deba@2382:         Node t = _graph.target(it);
deba@2382:         if (_comp->comp(s) == _comp->comp(t)) continue;
deba@2382: 
deba@2382:         Value cost = _dijkstra.dist(s) + _dijkstra.dist(t) + _cost[it];
deba@2382: 
deba@2382:         if (_comp->comp(s) < _comp->comp(t)) {
deba@2382:           externals.push_back(External(_comp->comp(s), _comp->comp(t), 
deba@2382:                                        it, cost));
deba@2382:         } else {
deba@2382:           externals.push_back(External(_comp->comp(t), _comp->comp(s), 
deba@2382:                                        it, cost));
deba@2382:         }
deba@2382:       }
deba@2382:       std::sort(externals.begin(), externals.end(), ExternalLess());
deba@2382: 
deba@2382:       SmartUGraph aux_graph;
deba@2382:       std::vector<SmartUGraph::Node> aux_nodes;
deba@2382: 
deba@2382:       for (int i = 0; i < _terminal_num; ++i) {
deba@2382:         aux_nodes.push_back(aux_graph.addNode());
deba@2382:       }
deba@2382: 
deba@2382:       SmartUGraph::UEdgeMap<Value> aux_cost(aux_graph);
deba@2382:       SmartUGraph::UEdgeMap<UEdge> cross(aux_graph);
deba@2382:       {
deba@2382:         int i = 0;
deba@2386:         while (i < int(externals.size())) {
deba@2382:           int sn = externals[i].source;
deba@2382:           int tn = externals[i].target;
deba@2382:           Value ev = externals[i].value;
deba@2382:           UEdge ee = externals[i].uedge;
deba@2382:           ++i;
deba@2386:           while (i < int(externals.size()) && 
deba@2382:                  sn == externals[i].source && tn == externals[i].target) {
deba@2382:             if (externals[i].value < ev) {
deba@2382:               ev = externals[i].value;
deba@2382:               ee = externals[i].uedge;
deba@2382:             }
deba@2382:             ++i;
deba@2382:           }
deba@2382:           SmartUGraph::UEdge ne = 
deba@2382:             aux_graph.addEdge(aux_nodes[sn], aux_nodes[tn]);
deba@2382:           aux_cost.set(ne, ev);
deba@2382:           cross.set(ne, ee);
deba@2382:         }
deba@2382:       }
deba@2382: 
deba@2382:       std::vector<SmartUGraph::UEdge> aux_tree_edges;
deba@2382:       BackInserterBoolMap<std::vector<SmartUGraph::UEdge> >
deba@2382:         aux_tree_map(aux_tree_edges);
deba@2382:       prim(aux_graph, aux_cost, aux_tree_map);
deba@2382: 
deba@2382:       for (std::vector<SmartUGraph::UEdge>::iterator 
deba@2382:              it = aux_tree_edges.begin(); it != aux_tree_edges.end(); ++it) {
deba@2382:         Node node;
deba@2382:         node = _graph.source(cross[*it]);
deba@2382:         while (node != INVALID && !(*_filter)[node]) {
deba@2382:           _filter->set(node, true);
deba@2382:           node = (*_pred)[node] != INVALID ? 
deba@2382:             _graph.source((*_pred)[node]) : INVALID;
deba@2382:         }
deba@2382:         node = _graph.target(cross[*it]);
deba@2382:         while (node != INVALID && !(*_filter)[node]) {
deba@2382:           _filter->set(node, true);
deba@2382:           node = (*_pred)[node] != INVALID ? 
deba@2382:             _graph.source((*_pred)[node]) : INVALID;
deba@2382:         }
deba@2382:       }
deba@2382: 
deba@2399:       _value = prim(nodeSubUGraphAdaptor(_graph, *_filter), _cost, *_tree);
deba@2382:             
deba@2382:     }
deba@2382: 
deba@2382:     /// \brief Checks if an edge is in the Steiner-tree or not.
deba@2382:     ///
deba@2382:     /// Checks if an edge is in the Steiner-tree or not.
deba@2382:     /// \param e is the edge that will be checked
deba@2382:     /// \return \c true if e is in the Steiner-tree, \c false otherwise
deba@2382:     bool tree(UEdge e){
deba@2382:       return (*_tree)[e];
deba@2382:     }
deba@2382: 
deba@2382:     /// \brief Checks if the node is in the Steiner-tree or not.
deba@2382:     ///
deba@2382:     /// Checks if a node is in the Steiner-tree or not.
deba@2382:     /// \param n is the node that will be checked
deba@2382:     /// \return \c true if n is in the Steiner-tree, \c false otherwise
deba@2382:     bool tree(Node n){
deba@2382:       return (*_filter)[n];
deba@2382:     }
deba@2399: 
deba@2399:     /// \brief Checks if the node is a Steiner-node.
deba@2399:     ///
deba@2399:     /// Checks if the node is a Steiner-node (i.e. a tree node but
deba@2399:     /// not terminal).
deba@2399:     /// \param n is the node that will be checked
deba@2399:     /// \return \c true if n is a Steiner-node, \c false otherwise
deba@2399:     bool steiner(Node n){
deba@2399:       return (*_filter)[n] && (*_pred)[n] != INVALID;
deba@2399:     }
deba@2399: 
deba@2399:     /// \brief Checks if the node is a terminal.
deba@2399:     ///
deba@2399:     /// Checks if the node is a terminal.
deba@2399:     /// \param n is the node that will be checked
deba@2399:     /// \return \c true if n is a terminal, \c false otherwise
deba@2399:     bool terminal(Node n){
deba@2399:       return _dijkstra.reached(n) && (*_pred)[n] == INVALID;
deba@2399:     }
deba@2382:     
deba@2399:     /// \brief The total cost of the tree
deba@2399:     ///
deba@2399:     /// The total cost of the constructed tree. The calculated value does
deba@2399:     /// not exceed the double of the optimal value.
deba@2399:     Value treeValue() const {
deba@2399:       return _value;
deba@2399:     }
deba@2399:     
deba@2382:   };
deba@2382: 
deba@2382: } //END OF NAMESPACE LEMON
deba@2382: 
deba@2382: #endif