78362ea3847d3bf9cb5a7caaa01199e7637ddd28
[libdai.git] / include / dai / weightedgraph.h
1 /* This file is part of libDAI - http://www.libdai.org/
2 *
3 * libDAI is licensed under the terms of the GNU General Public License version
4 * 2, or (at your option) any later version. libDAI is distributed without any
5 * warranty. See the file COPYING for more details.
6 *
7 * Copyright (C) 2006-2010 Joris Mooij [joris dot mooij at libdai dot org]
8 * Copyright (C) 2006-2007 Radboud University Nijmegen, The Netherlands
9 */
10
11
12 /** \file
13 * \brief Defines some utility functions for (weighted) undirected graphs, trees and rooted trees.
14 * \todo Improve general support for graphs and trees.
15 */
16
17
18 #ifndef __defined_libdai_weightedgraph_h
19 #define __defined_libdai_weightedgraph_h
20
21
22 #include <vector>
23 #include <map>
24 #include <iostream>
25 #include <set>
26 #include <limits>
27 #include <climits> // Work-around for bug in boost graph library
28 #include <dai/util.h>
29 #include <dai/exceptions.h>
30
31 #include <boost/graph/adjacency_list.hpp>
32 #include <boost/graph/prim_minimum_spanning_tree.hpp>
33 #include <boost/graph/kruskal_min_spanning_tree.hpp>
34
35
36 namespace dai {
37
38
39 /// Represents a directed edge
40 class DEdge {
41 public:
42 /// First node index (source of edge)
43 size_t first;
44
45 /// Second node index (target of edge)
46 size_t second;
47
48 /// Default constructor
49 DEdge() : first(0), second(0) {}
50
51 /// Constructs a directed edge pointing from \a m1 to \a m2
52 DEdge( size_t m1, size_t m2 ) : first(m1), second(m2) {}
53
54 /// Tests for equality
55 bool operator==( const DEdge &x ) const { return ((first == x.first) && (second == x.second)); }
56
57 /// Smaller-than operator (performs lexicographical comparison)
58 bool operator<( const DEdge &x ) const {
59 return( (first < x.first) || ((first == x.first) && (second < x.second)) );
60 }
61
62 /// Writes a directed edge to an output stream
63 friend std::ostream & operator << (std::ostream & os, const DEdge & e) {
64 os << "(" << e.first << "->" << e.second << ")";
65 return os;
66 }
67 };
68
69
70 /// Represents an undirected edge
71 class UEdge {
72 public:
73 /// First node index
74 size_t first;
75
76 /// Second node index
77 size_t second;
78
79 /// Default constructor
80 UEdge() : first(0), second(0) {}
81
82 /// Constructs an undirected edge between \a m1 and \a m2
83 UEdge( size_t m1, size_t m2 ) : first(m1), second(m2) {}
84
85 /// Construct from DEdge
86 UEdge( const DEdge &e ) : first(e.first), second(e.second) {}
87
88 /// Tests for inequality (disregarding the ordering of the nodes)
89 bool operator==( const UEdge &x ) {
90 return ((first == x.first) && (second == x.second)) || ((first == x.second) && (second == x.first));
91 }
92
93 /// Smaller-than operator
94 bool operator<( const UEdge &x ) const {
95 size_t s = std::min( first, second );
96 size_t l = std::max( first, second );
97 size_t xs = std::min( x.first, x.second );
98 size_t xl = std::max( x.first, x.second );
99 return( (s < xs) || ((s == xs) && (l < xl)) );
100 }
101
102 /// Writes an undirected edge to an output stream
103 friend std::ostream & operator << (std::ostream & os, const UEdge & e) {
104 if( e.first < e.second )
105 os << "{" << e.first << "--" << e.second << "}";
106 else
107 os << "{" << e.second << "--" << e.first << "}";
108 return os;
109 }
110 };
111
112
113 /// Represents an undirected graph, implemented as a std::set of undirected edges
114 class GraphEL : public std::set<UEdge> {
115 public:
116 /// Default constructor
117 GraphEL() {}
118
119 /// Construct from range of objects that can be cast to UEdge
120 template <class InputIterator>
121 GraphEL( InputIterator begin, InputIterator end ) {
122 insert( begin, end );
123 }
124 };
125
126
127 /// Represents an undirected weighted graph, with weights of type \a T, implemented as a std::map mapping undirected edges to weights
128 template<class T> class WeightedGraph : public std::map<UEdge, T> {};
129
130
131 /// Represents a rooted tree, implemented as a vector of directed edges
132 /** By convention, the edges are stored such that they point away from
133 * the root and such that edges nearer to the root come before edges
134 * farther away from the root.
135 */
136 class RootedTree : public std::vector<DEdge> {
137 public:
138 /// Default constructor
139 RootedTree() {}
140
141 /// Constructs a rooted tree from a tree and a root
142 /** \pre T has no cycles and contains node \a Root
143 */
144 RootedTree( const GraphEL &T, size_t Root );
145 };
146
147
148 /// Constructs a minimum spanning tree from the (non-negatively) weighted graph \a G.
149 /** \param G Weighted graph that should have non-negative weights.
150 * \param usePrim If true, use Prim's algorithm (complexity O(E log(V))), otherwise, use Kruskal's algorithm (complexity O(E log(E))).
151 * \note Uses implementation from Boost Graph Library.
152 * \note The vertices of \a G must be in the range [0,N) where N is the number of vertices of \a G.
153 */
154 template<typename T> RootedTree MinSpanningTree( const WeightedGraph<T> &G, bool usePrim ) {
155 RootedTree result;
156 if( G.size() > 0 ) {
157 using namespace boost;
158 using namespace std;
159 typedef adjacency_list< vecS, vecS, undirectedS, no_property, property<edge_weight_t, double> > boostGraph;
160
161 set<size_t> nodes;
162 vector<UEdge> edges;
163 vector<double> weights;
164 edges.reserve( G.size() );
165 weights.reserve( G.size() );
166 for( typename WeightedGraph<T>::const_iterator e = G.begin(); e != G.end(); e++ ) {
167 weights.push_back( e->second );
168 edges.push_back( e->first );
169 nodes.insert( e->first.first );
170 nodes.insert( e->first.second );
171 }
172
173 size_t N = nodes.size();
174 for( set<size_t>::const_iterator it = nodes.begin(); it != nodes.end(); it++ )
175 if( *it >= N )
176 DAI_THROWE(RUNTIME_ERROR,"Vertices must be in range [0..N) where N is the number of vertices.");
177
178 boostGraph g( edges.begin(), edges.end(), weights.begin(), nodes.size() );
179 size_t root = *(nodes.begin());
180 GraphEL tree;
181 if( usePrim ) {
182 // Prim's algorithm
183 vector< graph_traits< boostGraph >::vertex_descriptor > p(N);
184 prim_minimum_spanning_tree( g, &(p[0]) );
185
186 // Store tree edges in result
187 for( size_t i = 0; i != p.size(); i++ ) {
188 if( p[i] != i )
189 tree.insert( UEdge( p[i], i ) );
190 }
191 } else {
192 // Kruskal's algorithm
193 vector< graph_traits< boostGraph >::edge_descriptor > t;
194 t.reserve( N - 1 );
195 kruskal_minimum_spanning_tree( g, std::back_inserter(t) );
196
197 // Store tree edges in result
198 for( size_t i = 0; i != t.size(); i++ ) {
199 size_t v1 = source( t[i], g );
200 size_t v2 = target( t[i], g );
201 if( v1 != v2 )
202 tree.insert( UEdge( v1, v2 ) );
203 }
204 }
205
206 // Direct edges in order to obtain a rooted tree
207 result = RootedTree( tree, root );
208 }
209 return result;
210 }
211
212
213 /// Constructs a minimum spanning tree from the (non-negatively) weighted graph \a G.
214 /** \param G Weighted graph that should have non-negative weights.
215 * \param usePrim If true, use Prim's algorithm (complexity O(E log(V))), otherwise, use Kruskal's algorithm (complexity O(E log(E))).
216 * \note Uses implementation from Boost Graph Library.
217 * \note The vertices of \a G must be in the range [0,N) where N is the number of vertices of \a G.
218 */
219 template<typename T> RootedTree MaxSpanningTree( const WeightedGraph<T> &G, bool usePrim ) {
220 if( G.size() == 0 )
221 return RootedTree();
222 else {
223 T maxweight = G.begin()->second;
224 for( typename WeightedGraph<T>::const_iterator it = G.begin(); it != G.end(); it++ )
225 if( it->second > maxweight )
226 maxweight = it->second;
227 // make a copy of the graph
228 WeightedGraph<T> gr( G );
229 // invoke MinSpanningTree with negative weights
230 // (which have to be shifted to satisfy positivity criterion)
231 for( typename WeightedGraph<T>::iterator it = gr.begin(); it != gr.end(); it++ )
232 it->second = maxweight - it->second;
233 return MinSpanningTree( gr, usePrim );
234 }
235 }
236
237
238 } // end of namespace dai
239
240
241 #endif