[Benjamin Piwowarski] Renamed "foreach" macro into "bforeach" to avoid conflicts...
[libdai.git] / src / graph.cpp
1 /* This file is part of libDAI - http://www.libdai.org/
2 *
3 * Copyright (c) 2006-2011, The libDAI authors. All rights reserved.
4 *
5 * Use of this source code is governed by a BSD-style license that can be found in the LICENSE file.
6 */
7
8
9 #include <dai/graph.h>
10
11
12 namespace dai {
13
14
15 using namespace std;
16
17
18 GraphAL& GraphAL::addEdge( size_t n1, size_t n2, bool check ) {
19 DAI_ASSERT( n1 < nrNodes() );
20 DAI_ASSERT( n2 < nrNodes() );
21 bool exists = false;
22 if( check ) {
23 // Check whether the edge already exists
24 bforeach( const Neighbor &n, nb(n1) )
25 if( n == n2 ) {
26 exists = true;
27 break;
28 }
29 }
30 if( !exists && n1 != n2 ) { // Add edge
31 Neighbor nb_1( nb(n1).size(), n2, nb(n2).size() );
32 Neighbor nb_2( nb_1.dual, n1, nb_1.iter );
33 nb(n1).push_back( nb_1 );
34 nb(n2).push_back( nb_2 );
35 }
36 return *this;
37 }
38
39
40 void GraphAL::eraseNode( size_t n ) {
41 DAI_ASSERT( n < nrNodes() );
42 // Erase neighbor entry of node n
43 _nb.erase( _nb.begin() + n );
44 // Adjust neighbor entries of nodes
45 for( size_t n2 = 0; n2 < nrNodes(); n2++ ) {
46 for( size_t iter = 0; iter < nb(n2).size(); ) {
47 Neighbor &m = nb(n2, iter);
48 if( m.node == n ) {
49 // delete this entry, because it points to the deleted node
50 nb(n2).erase( nb(n2).begin() + iter );
51 } else {
52 // update this entry and the corresponding dual of the neighboring node
53 if( m.node > n )
54 m.node--;
55 nb( m.node, m.dual ).dual = iter;
56 m.iter = iter++;
57 }
58 }
59 }
60 }
61
62
63 void GraphAL::eraseEdge( size_t n1, size_t n2 ) {
64 DAI_ASSERT( n1 < nrNodes() );
65 DAI_ASSERT( n2 < nrNodes() );
66 size_t iter;
67 // Search for edge among neighbors of n1
68 for( iter = 0; iter < nb(n1).size(); iter++ )
69 if( nb(n1, iter).node == n2 ) {
70 // Remove it
71 nb(n1).erase( nb(n1).begin() + iter );
72 break;
73 }
74 // Change the iter and dual values of the subsequent neighbors
75 for( ; iter < nb(n1).size(); iter++ ) {
76 Neighbor &m = nb( n1, iter );
77 m.iter = iter;
78 nb( m.node, m.dual ).dual = iter;
79 }
80 // Search for edge among neighbors of n2
81 for( iter = 0; iter < nb(n2).size(); iter++ )
82 if( nb(n2, iter).node == n1 ) {
83 // Remove it
84 nb(n2).erase( nb(n2).begin() + iter );
85 break;
86 }
87 // Change the iter and node values of the subsequent neighbors
88 for( ; iter < nb(n2).size(); iter++ ) {
89 Neighbor &m = nb( n2, iter );
90 m.iter = iter;
91 nb( m.node, m.dual ).dual = iter;
92 }
93 }
94
95
96 SmallSet<size_t> GraphAL::nbSet( size_t n ) const {
97 SmallSet<size_t> result;
98 bforeach( const Neighbor &m, nb(n) )
99 result |= m;
100 return result;
101 }
102
103
104 bool GraphAL::isConnected() const {
105 if( nrNodes() == 0 ) {
106 return true;
107 } else {
108 std::vector<bool> incomponent( nrNodes(), false );
109
110 incomponent[0] = true;
111 bool found_new_nodes;
112 do {
113 found_new_nodes = false;
114
115 // For all nodes, check if they are connected with the (growing) component
116 for( size_t n1 = 0; n1 < nrNodes(); n1++ )
117 if( !incomponent[n1] ) {
118 bforeach( const Neighbor &n2, nb(n1) ) {
119 if( incomponent[n2] ) {
120 found_new_nodes = true;
121 incomponent[n1] = true;
122 break;
123 }
124 }
125 }
126 } while( found_new_nodes );
127
128 // Check if there are remaining nodes (not in the component)
129 bool all_connected = true;
130 for( size_t n1 = 0; (n1 < nrNodes()) && all_connected; n1++ )
131 if( !incomponent[n1] )
132 all_connected = false;
133
134 return all_connected;
135
136 // BGL implementation is slower...
137 /* using namespace boost;
138 typedef adjacency_list< vecS, vecS, undirectedS, property<vertex_distance_t, int> > boostGraphAL;
139 typedef pair<size_t, size_t> E;
140
141 // Copy graph structure into boostGraphAL object
142 vector<E> edges;
143 edges.reserve( nrEdges() );
144 for( size_t n1 = 0; n1 < nrNodes(); n1++ )
145 bforeach( const Neighbor &n2, nb(n1) )
146 if( n1 < n2 )
147 edges.push_back( E( n1, n2 ) );
148 boostGraphAL g( edges.begin(), edges.end(), nrNodes() );
149
150 // Construct connected components using Boost GraphAL Library
151 std::vector<int> component( num_vertices( g ) );
152 int num_comp = connected_components( g, make_iterator_property_map(component.begin(), get(vertex_index, g)) );
153
154 return (num_comp == 1);
155 */
156 }
157 }
158
159
160 bool GraphAL::isTree() const {
161 typedef vector<Edge> levelType; // first is node, second is its parent
162 vector<levelType> levels;
163
164 if( nrNodes() == 0 )
165 return true;
166 else {
167 // start with root node 0
168 levels.push_back( levelType( 1, Edge( 0, 0 ) ) );
169 size_t treeSize = 1;
170 bool foundCycle = false;
171 do {
172 levels.push_back( levelType() );
173 const levelType &prevLevel = levels[levels.size() - 2];
174 // build new level: add all neighbors of nodes in the previous level
175 // (without backtracking), aborting if a cycle is detected
176 for( size_t e = 0; e < prevLevel.size(); e++ ) {
177 size_t n2 = prevLevel[e].first; // for all nodes n2 in the previous level
178 bforeach( const Neighbor &n1, nb(n2) ) { // for all neighbors n1 of n2
179 if( n1 != prevLevel[e].second ) { // no backtracking allowed
180 for( size_t l = 0; l < levels.size() && !foundCycle; l++ )
181 for( size_t f = 0; f < levels[l].size() && !foundCycle; f++ )
182 if( levels[l][f].first == n1 )
183 // n1 has been visited before -> found a cycle
184 foundCycle = true;
185 if( !foundCycle )
186 // add n1 (and its parent n2) to current level
187 levels.back().push_back( Edge( n1, n2 ) );
188 }
189 if( foundCycle )
190 break;
191 }
192 if( foundCycle )
193 break;
194 }
195 treeSize += levels.back().size();
196 } while( (levels.back().size() != 0) && !foundCycle );
197 if( treeSize == nrNodes() && !foundCycle )
198 return true;
199 else
200 return false;
201 }
202 }
203
204
205 void GraphAL::printDot( std::ostream& os ) const {
206 os << "graph GraphAL {" << endl;
207 os << "node[shape=circle,width=0.4,fixedsize=true];" << endl;
208 for( size_t n = 0; n < nrNodes(); n++ )
209 os << "\tx" << n << ";" << endl;
210 for( size_t n1 = 0; n1 < nrNodes(); n1++ )
211 bforeach( const Neighbor &n2, nb(n1) )
212 if( n1 < n2 )
213 os << "\tx" << n1 << " -- x" << n2 << ";" << endl;
214 os << "}" << endl;
215 }
216
217
218 void GraphAL::checkConsistency() const {
219 size_t N = nrNodes();
220 for( size_t n1 = 0; n1 < N; n1++ ) {
221 size_t iter = 0;
222 bforeach( const Neighbor &n2, nb(n1) ) {
223 DAI_ASSERT( n2.iter == iter );
224 DAI_ASSERT( n2.node < N );
225 DAI_ASSERT( n2.dual < nb(n2).size() );
226 DAI_ASSERT( nb(n2, n2.dual) == n1 );
227 iter++;
228 }
229 }
230 }
231
232
233 GraphAL createGraphFull( size_t N ) {
234 GraphAL result( N );
235 for( size_t i = 0; i < N; i++ )
236 for( size_t j = i+1; j < N; j++ )
237 result.addEdge( i, j, false );
238 return result;
239 }
240
241
242 GraphAL createGraphGrid( size_t N1, size_t N2, bool periodic ) {
243 GraphAL result( N1*N2 );
244 if( N1 == 1 && N2 == 1 )
245 return result;
246 for( size_t i1 = 0; i1 < N1; i1++ )
247 for( size_t i2 = 0; i2 < N2; i2++ ) {
248 if( i1+1 < N1 || periodic )
249 result.addEdge( i1*N2 + i2, ((i1+1)%N1)*N2 + i2, N1 <= 2 );
250 if( i2+1 < N2 || periodic )
251 result.addEdge( i1*N2 + i2, i1*N2 + ((i2+1)%N2), N2 <= 2 );
252 }
253 return result;
254 }
255
256
257 GraphAL createGraphGrid3D( size_t N1, size_t N2, size_t N3, bool periodic ) {
258 GraphAL result( N1*N2*N3 );
259 for( size_t i1 = 0; i1 < N1; i1++ )
260 for( size_t i2 = 0; i2 < N2; i2++ )
261 for( size_t i3 = 0; i3 < N3; i3++ ) {
262 if( i1+1 < N1 || periodic )
263 result.addEdge( i1*N2*N3 + i2*N3 + i3, ((i1+1)%N1)*N2*N3 + i2*N3 + i3, N1 <= 2 );
264 if( i2+1 < N2 || periodic )
265 result.addEdge( i1*N2*N3 + i2*N3 + i3, i1*N2*N3 + ((i2+1)%N2)*N3 + i3, N2 <= 2 );
266 if( i3+1 < N3 || periodic )
267 result.addEdge( i1*N2*N3 + i2*N3 + i3, i1*N2*N3 + i2*N3 + ((i3+1)%N3), N3 <= 2 );
268 }
269 return result;
270 }
271
272
273 GraphAL createGraphLoop( size_t N ) {
274 GraphAL result( N );
275 for( size_t i = 0; i < N; i++ )
276 result.addEdge( i, (i+1)%N, N <= 2 );
277 return result;
278 }
279
280
281 GraphAL createGraphTree( size_t N ) {
282 GraphAL result( N );
283 for( size_t i = 1; i < N; i++ ) {
284 size_t j = rnd_int( 0, i-1 );
285 result.addEdge( i, j, false );
286 }
287 return result;
288 }
289
290
291 GraphAL createGraphRegular( size_t N, size_t d ) {
292 DAI_ASSERT( (N * d) % 2 == 0 );
293 DAI_ASSERT( d < N );
294
295 GraphAL G( N );
296 if( d > 0 ) {
297 bool ready = false;
298 size_t tries = 0;
299 while( !ready ) {
300 tries++;
301
302 // Start with N*d points {0,1,...,N*d-1} (N*d even) in N groups.
303 // Put U = {0,1,...,N*d-1}. (U denotes the set of unpaired points.)
304 vector<size_t> U;
305 U.reserve( N * d );
306 for( size_t i = 0; i < N * d; i++ )
307 U.push_back( i );
308
309 // Repeat the following until no suitable pair can be found: Choose
310 // two random points i and j in U, and if they are suitable, pair
311 // i with j and delete i and j from U.
312 G = GraphAL( N );
313 bool finished = false;
314 while( !finished ) {
315 random_shuffle( U.begin(), U.end(), rnd );
316 size_t i1, i2;
317 bool suit_pair_found = false;
318 for( i1 = 0; i1 < U.size()-1 && !suit_pair_found; i1++ )
319 for( i2 = i1+1; i2 < U.size() && !suit_pair_found; i2++ )
320 if( ((U[i1] / d) != (U[i2] / d)) && !G.hasEdge( U[i1] / d, U[i2] / d ) ) {
321 // they are suitable
322 suit_pair_found = true;
323 G.addEdge( U[i1] / d, U[i2] / d, false );
324 U.erase( U.begin() + i2 ); // first remove largest
325 U.erase( U.begin() + i1 ); // remove smallest
326 }
327 if( !suit_pair_found || U.empty() )
328 finished = true;
329 }
330
331 if( U.empty() ) {
332 // G is a graph with edge from vertex r to vertex s if and only if
333 // there is a pair containing points in the r'th and s'th groups.
334 // If G is d-regular, output, otherwise return to Step 1.
335 ready = true;
336 for( size_t n = 0; n < N; n++ )
337 if( G.nb(n).size() != d ) {
338 ready = false;
339 break;
340 }
341 } else
342 ready = false;
343 }
344 }
345
346 return G;
347 }
348
349
350 } // end of namespace dai