Merge branch 'vaske'
[libdai.git] / include / dai / weightedgraph.h
1 /* Copyright (C) 2006-2008 Joris Mooij [joris dot mooij at tuebingen dot mpg dot de]
2 Radboud University Nijmegen, The Netherlands /
3 Max Planck Institute for Biological Cybernetics, Germany
4
5 This file is part of libDAI.
6
7 libDAI is free software; you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation; either version 2 of the License, or
10 (at your option) any later version.
11
12 libDAI is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with libDAI; if not, write to the Free Software
19 Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
20 */
21
22
23 /** \file
24 * \brief Defines some utility functions for weighted graphs
25 * \todo Improve documentation
26 * \todo Improve general support for graphs and trees.
27 */
28
29
30 #ifndef __defined_libdai_weightedgraph_h
31 #define __defined_libdai_weightedgraph_h
32
33
34 #include <vector>
35 #include <map>
36 #include <iostream>
37 #include <set>
38 #include <cassert>
39 #include <limits>
40 #include <climits> // Work-around for bug in boost graph library
41
42 #include <boost/graph/adjacency_list.hpp>
43 #include <boost/graph/prim_minimum_spanning_tree.hpp>
44
45
46 namespace dai {
47
48
49 /// Represents a directed edge pointing from n1 to n2
50 class DEdge {
51 public:
52 size_t n1; ///< First node index
53 size_t n2; ///< Second node index
54
55 /// Default constructor
56 DEdge() {}
57
58 /// Constructor
59 DEdge( size_t m1, size_t m2 ) : n1(m1), n2(m2) {}
60
61 /// Tests for equality
62 bool operator==( const DEdge &x ) const { return ((n1 == x.n1) && (n2 == x.n2)); }
63
64 /// Tests for inequality
65 bool operator!=( const DEdge &x ) const { return !(*this == x); }
66
67 /// Smaller-than operator (performs lexicographical comparison)
68 bool operator<( const DEdge &x ) const {
69 return( (n1 < x.n1) || ((n1 == x.n1) && (n2 < x.n2)) );
70 }
71
72 /// Writes a DEdge to an output stream
73 friend std::ostream & operator << (std::ostream & os, const DEdge & e) {
74 os << "(" << e.n1 << "," << e.n2 << ")";
75 return os;
76 }
77 };
78
79
80 /// Undirected edge between nodes n1 and n2
81 class UEdge {
82 public:
83 size_t n1; ///< First node index
84 size_t n2; ///< Second node index
85
86 /// Default constructor
87 UEdge() {}
88
89 /// Constructor
90 UEdge( size_t m1, size_t m2 ) : n1(m1), n2(m2) {}
91
92 /// Construct from DEdge
93 UEdge( const DEdge & e ) : n1(e.n1), n2(e.n2) {}
94
95 /// Tests for inequality (disregarding the ordering of n1 and n2)
96 bool operator==( const UEdge &x ) {
97 return ((n1 == x.n1) && (n2 == x.n2)) || ((n1 == x.n2) && (n2 == x.n1));
98 }
99
100 /// Smaller-than operator
101 bool operator<( const UEdge &x ) const {
102 size_t s = n1, l = n2;
103 if( s > l )
104 std::swap( s, l );
105 size_t xs = x.n1, xl = x.n2;
106 if( xs > xl )
107 std::swap( xs, xl );
108 return( (s < xs) || ((s == xs) && (l < xl)) );
109 }
110
111 /// Writes a UEdge to an output stream
112 friend std::ostream & operator << (std::ostream & os, const UEdge & e) {
113 if( e.n1 < e.n2 )
114 os << "{" << e.n1 << "," << e.n2 << "}";
115 else
116 os << "{" << e.n2 << "," << e.n1 << "}";
117 return os;
118 }
119 };
120
121
122 /// Vector of UEdge
123 typedef std::vector<UEdge> UEdgeVec;
124
125 /// Vector of DEdge
126 typedef std::vector<DEdge> DEdgeVec;
127
128 /// Represents an undirected weighted graph, with weights of type T
129 template<class T> class WeightedGraph : public std::map<UEdge, T> {};
130
131 /// Represents an undirected graph
132 typedef std::set<UEdge> Graph;
133
134
135 /// Uses Prim's algorithm to construct a minimal spanning tree from the (positively) weighted graph G.
136 /** Uses implementation in Boost Graph Library.
137 */
138 template<typename T> DEdgeVec MinSpanningTreePrims( const WeightedGraph<T> &G ) {
139 DEdgeVec result;
140 if( G.size() > 0 ) {
141 using namespace boost;
142 using namespace std;
143 typedef adjacency_list< vecS, vecS, undirectedS, property<vertex_distance_t, int>, property<edge_weight_t, double> > boostGraph;
144 typedef pair<size_t, size_t> E;
145
146 set<size_t> nodes;
147 vector<E> edges;
148 vector<double> weights;
149 edges.reserve( G.size() );
150 weights.reserve( G.size() );
151 for( typename WeightedGraph<T>::const_iterator e = G.begin(); e != G.end(); e++ ) {
152 weights.push_back( e->second );
153 edges.push_back( E( e->first.n1, e->first.n2 ) );
154 nodes.insert( e->first.n1 );
155 nodes.insert( e->first.n2 );
156 }
157
158 boostGraph g( edges.begin(), edges.end(), weights.begin(), nodes.size() );
159 vector< graph_traits< boostGraph >::vertex_descriptor > p( num_vertices(g) );
160 prim_minimum_spanning_tree( g, &(p[0]) );
161
162 // Store tree edges in result
163 result.reserve( nodes.size() - 1 );
164 size_t root = 0;
165 for( size_t i = 0; i != p.size(); i++ )
166 if( p[i] != i )
167 result.push_back( DEdge( p[i], i ) );
168 else
169 root = i;
170
171 // We have to store the minimum spanning tree in the right
172 // order, such that for all (i1, j1), (i2, j2) in result,
173 // if j1 == i2 then (i1, j1) comes before (i2, j2) in result.
174 // We do this by reordering the contents of result, effectively
175 // growing the tree starting at the root. At each step,
176 // result[0..N-1] are the edges already added to the tree,
177 // whereas the other elements of result still have to be added.
178 // The elements of nodes are the vertices that still have to
179 // be added to the tree.
180
181 // Start with the root
182 nodes.erase( root );
183 size_t N = 0;
184
185 // Iteratively add edges and nodes to the growing tree
186 while( N != result.size() ) {
187 for( size_t e = N; e != result.size(); e++ ) {
188 bool e1_in_tree = !nodes.count( result[e].n1 );
189 if( e1_in_tree ) {
190 nodes.erase( result[e].n2 );
191 swap( result[N], result[e] );
192 N++;
193 break;
194 }
195 }
196 }
197 }
198
199 return result;
200 }
201
202
203 /// Use Prim's algorithm to construct a minimal spanning tree from the (positively) weighted graph G.
204 /** Uses implementation in Boost Graph Library.
205 */
206 template<typename T> DEdgeVec MaxSpanningTreePrims( const WeightedGraph<T> & Graph ) {
207 if( Graph.size() == 0 )
208 return DEdgeVec();
209 else {
210 T maxweight = Graph.begin()->second;
211 for( typename WeightedGraph<T>::const_iterator it = Graph.begin(); it != Graph.end(); it++ )
212 if( it->second > maxweight )
213 maxweight = it->second;
214 // make a copy of the graph
215 WeightedGraph<T> gr( Graph );
216 // invoke MinSpanningTreePrims with negative weights
217 // (which have to be shifted to satisfy positivity criterion)
218 for( typename WeightedGraph<T>::iterator it = gr.begin(); it != gr.end(); it++ )
219 it->second = maxweight - it->second;
220 return MinSpanningTreePrims( gr );
221 }
222 }
223
224
225 /// Constructs a rooted tree from a tree and a root
226 DEdgeVec GrowRootedTree( const Graph & T, size_t Root );
227
228
229 /// Constructs a random undirected graph of N nodes, where each node has connectivity d
230 UEdgeVec RandomDRegularGraph( size_t N, size_t d );
231
232
233 } // end of namespace dai
234
235
236 #endif