Added DAG class and various minor improvements
[libdai.git] / src / bipgraph.cpp
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 #include <dai/bipgraph.h>
13
14
15 namespace dai {
16
17
18 using namespace std;
19
20
21 BipartiteGraph& BipartiteGraph::addEdge( size_t n1, size_t n2, bool check ) {
22 DAI_ASSERT( n1 < nrNodes1() );
23 DAI_ASSERT( n2 < nrNodes2() );
24 bool exists = false;
25 if( check ) {
26 // Check whether the edge already exists
27 foreach( const Neighbor &nb2, nb1(n1) )
28 if( nb2 == n2 ) {
29 exists = true;
30 break;
31 }
32 }
33 if( !exists ) { // Add edge
34 Neighbor nb_1( nb1(n1).size(), n2, nb2(n2).size() );
35 Neighbor nb_2( nb_1.dual, n1, nb_1.iter );
36 nb1(n1).push_back( nb_1 );
37 nb2(n2).push_back( nb_2 );
38 }
39 return *this;
40 }
41
42
43 void BipartiteGraph::eraseNode1( size_t n1 ) {
44 DAI_ASSERT( n1 < nrNodes1() );
45 // Erase neighbor entry of node n1
46 _nb1.erase( _nb1.begin() + n1 );
47 // Adjust neighbor entries of nodes of type 2
48 for( size_t n2 = 0; n2 < nrNodes2(); n2++ ) {
49 for( size_t iter = 0; iter < nb2(n2).size(); ) {
50 Neighbor &m1 = nb2(n2, iter);
51 if( m1.node == n1 ) {
52 // delete this entry, because it points to the deleted node
53 nb2(n2).erase( nb2(n2).begin() + iter );
54 } else {
55 // update this entry and the corresponding dual of the neighboring node of type 1
56 if( m1.node > n1 )
57 m1.node--;
58 nb1( m1.node, m1.dual ).dual = iter;
59 m1.iter = iter++;
60 }
61 }
62 }
63 }
64
65
66 void BipartiteGraph::eraseNode2( size_t n2 ) {
67 DAI_ASSERT( n2 < nrNodes2() );
68 // Erase neighbor entry of node n2
69 _nb2.erase( _nb2.begin() + n2 );
70 // Adjust neighbor entries of nodes of type 1
71 for( size_t n1 = 0; n1 < nrNodes1(); n1++ ) {
72 for( size_t iter = 0; iter < nb1(n1).size(); ) {
73 Neighbor &m2 = nb1(n1, iter);
74 if( m2.node == n2 ) {
75 // delete this entry, because it points to the deleted node
76 nb1(n1).erase( nb1(n1).begin() + iter );
77 } else {
78 // update this entry and the corresponding dual of the neighboring node of type 2
79 if( m2.node > n2 )
80 m2.node--;
81 nb2( m2.node, m2.dual ).dual = iter;
82 m2.iter = iter++;
83 }
84 }
85 }
86 }
87
88
89 void BipartiteGraph::eraseEdge( size_t n1, size_t n2 ) {
90 DAI_ASSERT( n1 < nrNodes1() );
91 DAI_ASSERT( n2 < nrNodes2() );
92 size_t iter;
93 // Search for edge among neighbors of n1
94 for( iter = 0; iter < nb1(n1).size(); iter++ )
95 if( nb1(n1, iter).node == n2 ) {
96 // Remove it
97 nb1(n1).erase( nb1(n1).begin() + iter );
98 break;
99 }
100 // Change the iter and dual values of the subsequent neighbors
101 for( ; iter < nb1(n1).size(); iter++ ) {
102 Neighbor &m2 = nb1( n1, iter );
103 m2.iter = iter;
104 nb2( m2.node, m2.dual ).dual = iter;
105 }
106 // Search for edge among neighbors of n2
107 for( iter = 0; iter < nb2(n2).size(); iter++ )
108 if( nb2(n2, iter).node == n1 ) {
109 // Remove it
110 nb2(n2).erase( nb2(n2).begin() + iter );
111 break;
112 }
113 // Change the iter and node values of the subsequent neighbors
114 for( ; iter < nb2(n2).size(); iter++ ) {
115 Neighbor &m1 = nb2( n2, iter );
116 m1.iter = iter;
117 nb1( m1.node, m1.dual ).dual = iter;
118 }
119 }
120
121
122 SmallSet<size_t> BipartiteGraph::nb1Set( size_t n1 ) const {
123 SmallSet<size_t> result;
124 foreach( const Neighbor &n2, nb1(n1) )
125 result |= n2;
126 return result;
127 }
128
129
130 SmallSet<size_t> BipartiteGraph::nb2Set( size_t n2 ) const {
131 SmallSet<size_t> result;
132 foreach( const Neighbor &n1, nb2(n2) )
133 result |= n1;
134 return result;
135 }
136
137
138 SmallSet<size_t> BipartiteGraph::delta1( size_t n1, bool include ) const {
139 // get all second-order neighbors
140 SmallSet<size_t> result;
141 foreach( const Neighbor &n2, nb1(n1) )
142 foreach( const Neighbor &m1, nb2(n2) )
143 if( include || (m1 != n1) )
144 result |= m1;
145 return result;
146 }
147
148
149 SmallSet<size_t> BipartiteGraph::delta2( size_t n2, bool include ) const {
150 // store all second-order neighbors
151 SmallSet<size_t> result;
152 foreach( const Neighbor &n1, nb2(n2) )
153 foreach( const Neighbor &m2, nb1(n1) )
154 if( include || (m2 != n2) )
155 result |= m2;
156 return result;
157 }
158
159
160 bool BipartiteGraph::isConnected() const {
161 if( nrNodes1() == 0 ) {
162 return true;
163 } else {
164 /*
165 // The BGL implementation is significantly slower...
166 using namespace boost;
167 typedef adjacency_list< vecS, vecS, undirectedS, property<vertex_distance_t, int> > boostGraph;
168 typedef pair<size_t, size_t> E;
169
170 // Copy graph structure into boostGraph object
171 size_t N = nrNodes1();
172 vector<E> edges;
173 edges.reserve( nrEdges() );
174 for( size_t n1 = 0; n1 < nrNodes1(); n1++ )
175 foreach( const Neighbor &n2, nb1(n1) )
176 edges.push_back( E( n1, n2.node + N ) );
177 boostGraph g( edges.begin(), edges.end(), nrNodes1() + nrNodes2() );
178
179 // Construct connected components using Boost Graph Library
180 std::vector<int> component( num_vertices( g ) );
181 int num_comp = connected_components( g, make_iterator_property_map(component.begin(), get(vertex_index, g)) );
182
183 return (num_comp == 1);
184 */
185
186 std::vector<bool> incomponent1( nrNodes1(), false );
187 std::vector<bool> incomponent2( nrNodes2(), false );
188
189 incomponent1[0] = true;
190 bool found_new_nodes;
191 do {
192 found_new_nodes = false;
193
194 // For all nodes of type 2, check if they are connected with the (growing) component
195 for( size_t n2 = 0; n2 < nrNodes2(); n2++ )
196 if( !incomponent2[n2] ) {
197 foreach( const Neighbor &n1, nb2(n2) ) {
198 if( incomponent1[n1] ) {
199 found_new_nodes = true;
200 incomponent2[n2] = true;
201 break;
202 }
203 }
204 }
205
206 // For all nodes of type 1, check if they are connected with the (growing) component
207 for( size_t n1 = 0; n1 < nrNodes1(); n1++ )
208 if( !incomponent1[n1] ) {
209 foreach( const Neighbor &n2, nb1(n1) ) {
210 if( incomponent2[n2] ) {
211 found_new_nodes = true;
212 incomponent1[n1] = true;
213 break;
214 }
215 }
216 }
217 } while( found_new_nodes );
218
219 // Check if there are remaining nodes (not in the component)
220 bool all_connected = true;
221 for( size_t n1 = 0; (n1 < nrNodes1()) && all_connected; n1++ )
222 if( !incomponent1[n1] )
223 all_connected = false;
224 for( size_t n2 = 0; (n2 < nrNodes2()) && all_connected; n2++ )
225 if( !incomponent2[n2] )
226 all_connected = false;
227
228 return all_connected;
229 }
230 }
231
232
233 bool BipartiteGraph::isTree() const {
234 vector<levelType> levels;
235
236 bool foundCycle = false;
237 size_t nr_1 = 0;
238 size_t nr_2 = 0;
239
240 if( nrNodes1() == 0 )
241 return (nrNodes2() < 2 );
242 else if( nrNodes2() == 0 )
243 return (nrNodes1() < 2 );
244 else {
245 levelType newLevel;
246 do {
247 newLevel.ind1.clear();
248 newLevel.ind2.clear();
249 if( levels.size() == 0 ) {
250 size_t n1 = 0;
251 // add n1 to ind1
252 newLevel.ind1 = vector<size_t>( 1, n1 );
253 // add all neighbors of n1 to ind2
254 newLevel.ind2.reserve( nb1(n1).size() );
255 foreach( const Neighbor &n2, nb1(n1) )
256 newLevel.ind2.push_back( n2 );
257 } else {
258 const levelType &prevLevel = levels.back();
259 // build newLevel.ind1
260 foreach( size_t n2, prevLevel.ind2 ) { // for all n2 in the previous level
261 foreach( const Neighbor &n1, nb2(n2) ) { // for all neighbors n1 of n2
262 if( find( prevLevel.ind1.begin(), prevLevel.ind1.end(), n1 ) == prevLevel.ind1.end() ) { // n1 not in previous level
263 if( find( newLevel.ind1.begin(), newLevel.ind1.end(), n1 ) != newLevel.ind1.end() )
264 foundCycle = true; // n1 already in new level: we found a cycle
265 else
266 newLevel.ind1.push_back( n1 ); // add n1 to new level
267 }
268 if( foundCycle )
269 break;
270 }
271 if( foundCycle )
272 break;
273 }
274 // build newLevel.ind2
275 foreach( size_t n1, newLevel.ind1 ) { // for all n1 in this level
276 foreach( const Neighbor &n2, nb1(n1) ) { // for all neighbors n2 of n1
277 if( find( prevLevel.ind2.begin(), prevLevel.ind2.end(), n2 ) == prevLevel.ind2.end() ) { // n2 not in previous level
278 if( find( newLevel.ind2.begin(), newLevel.ind2.end(), n2 ) != newLevel.ind2.end() )
279 foundCycle = true; // n2 already in new level: we found a cycle
280 else
281 newLevel.ind2.push_back( n2 ); // add n2 to new level
282 }
283 if( foundCycle )
284 break;
285 }
286 if( foundCycle )
287 break;
288 }
289 }
290 levels.push_back( newLevel );
291 nr_1 += newLevel.ind1.size();
292 nr_2 += newLevel.ind2.size();
293 } while( ((newLevel.ind1.size() != 0) || (newLevel.ind2.size() != 0)) && !foundCycle );
294 if( nr_1 == nrNodes1() && nr_2 == nrNodes2() && !foundCycle )
295 return true;
296 else
297 return false;
298 }
299 }
300
301
302 void BipartiteGraph::printDot( std::ostream& os ) const {
303 os << "graph BipartiteGraph {" << endl;
304 os << "node[shape=circle,width=0.4,fixedsize=true];" << endl;
305 for( size_t n1 = 0; n1 < nrNodes1(); n1++ )
306 os << "\tx" << n1 << ";" << endl;
307 os << "node[shape=box,width=0.3,height=0.3,fixedsize=true];" << endl;
308 for( size_t n2 = 0; n2 < nrNodes2(); n2++ )
309 os << "\ty" << n2 << ";" << endl;
310 for( size_t n1 = 0; n1 < nrNodes1(); n1++ )
311 foreach( const Neighbor &n2, nb1(n1) )
312 os << "\tx" << n1 << " -- y" << n2 << ";" << endl;
313 os << "}" << endl;
314 }
315
316
317 void BipartiteGraph::checkConsistency() const {
318 size_t N1 = nrNodes1();
319 size_t N2 = nrNodes2();
320 for( size_t n1 = 0; n1 < N1; n1++ ) {
321 size_t iter = 0;
322 foreach( const Neighbor &n2, nb1(n1) ) {
323 DAI_ASSERT( n2.iter == iter );
324 DAI_ASSERT( n2.node < N2 );
325 DAI_ASSERT( n2.dual < nb2(n2).size() );
326 DAI_ASSERT( nb2(n2, n2.dual) == n1 );
327 iter++;
328 }
329 }
330 for( size_t n2 = 0; n2 < N2; n2++ ) {
331 size_t iter = 0;
332 foreach( const Neighbor &n1, nb2(n2) ) {
333 DAI_ASSERT( n1.iter == iter );
334 DAI_ASSERT( n1.node < N1 );
335 DAI_ASSERT( n1.dual < nb1(n1).size() );
336 DAI_ASSERT( nb1(n1, n1.dual) == n2 );
337 iter++;
338 }
339 }
340 }
341
342
343 } // end of namespace dai