1 /* -*- mode: C -*- */
2 /*
3 IGraph library.
4 Copyright (C) 2008-2012 Gabor Csardi <csardi.gabor@gmail.com>
5 334 Harvard street, Cambridge, MA 02139 USA
6
7 This program 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 This program 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 this program; if not, write to the Free Software
19 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
20 02110-1301 USA
21
22 */
23
24 #include "igraph_structural.h"
25 #include "igraph_error.h"
26 #include "igraph_adjlist.h"
27 #include "igraph_interface.h"
28
29 /**
30 * \function igraph_maximum_cardinality_search
31 * Maximum cardinality search
32 *
33 * This function implements the maximum cardinality search algorithm
34 * discussed in
35 * Robert E Tarjan and Mihalis Yannakakis: Simple linear-time
36 * algorithms to test chordality of graphs, test acyclicity of
37 * hypergraphs, and selectively reduce acyclic hypergraphs.
38 * SIAM Journal of Computation 13, 566--579, 1984.
39 *
40 * \param graph The input graph, which should be undirected and simple.
41 * of the edges is ignored.
42 * \param alpha Pointer to an initialized vector, the result is stored here.
43 * It will be resized, as needed. Upon return it contains
44 * the rank of the each vertex.
45 * \param alpham1 Pointer to an initialized vector or a \c NULL
46 * pointer. If not \c NULL, then the inverse of \p alpha is stored
47 * here.
48 * \return Error code.
49 *
50 * Time complexity: O(|V|+|E|), linear in terms of the number of
51 * vertices and edges.
52 *
53 * \sa \ref igraph_is_chordal().
54 */
55
igraph_maximum_cardinality_search(const igraph_t * graph,igraph_vector_t * alpha,igraph_vector_t * alpham1)56 int igraph_maximum_cardinality_search(const igraph_t *graph,
57 igraph_vector_t *alpha,
58 igraph_vector_t *alpham1) {
59
60 long int no_of_nodes = igraph_vcount(graph);
61 igraph_vector_long_t size;
62 igraph_vector_long_t head, next, prev; /* doubly linked list with head */
63 long int i;
64 igraph_adjlist_t adjlist;
65 igraph_bool_t simple;
66
67 /***************/
68 /* local j, v; */
69 /***************/
70
71 long int j, v;
72
73 if (igraph_is_directed(graph)) {
74 IGRAPH_ERROR("Maximum cardinality search works on undirected graphs only", IGRAPH_EINVAL);
75 }
76
77 igraph_is_simple(graph, &simple);
78 if (!simple) {
79 IGRAPH_ERROR("Maximum cardinality search works on simple graphs only", IGRAPH_EINVAL);
80 }
81
82 if (no_of_nodes == 0) {
83 igraph_vector_clear(alpha);
84 if (alpham1) {
85 igraph_vector_clear(alpham1);
86 }
87 return IGRAPH_SUCCESS;
88 }
89
90 IGRAPH_CHECK(igraph_vector_long_init(&size, no_of_nodes));
91 IGRAPH_FINALLY(igraph_vector_long_destroy, &size);
92 IGRAPH_CHECK(igraph_vector_long_init(&head, no_of_nodes));
93 IGRAPH_FINALLY(igraph_vector_long_destroy, &head);
94 IGRAPH_CHECK(igraph_vector_long_init(&next, no_of_nodes));
95 IGRAPH_FINALLY(igraph_vector_long_destroy, &next);
96 IGRAPH_CHECK(igraph_vector_long_init(&prev, no_of_nodes));
97 IGRAPH_FINALLY(igraph_vector_long_destroy, &prev);
98
99 IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_ALL));
100 IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
101
102 IGRAPH_CHECK(igraph_vector_resize(alpha, no_of_nodes));
103 if (alpham1) {
104 IGRAPH_CHECK(igraph_vector_resize(alpham1, no_of_nodes));
105 }
106
107 /***********************************************/
108 /* for i in [0,n-1] -> set(i) := emptyset rof; */
109 /***********************************************/
110
111 /* nothing to do, 'head' contains all zeros */
112
113 /*********************************************************/
114 /* for v in vertices -> size(v):=0; add v to set(0) rof; */
115 /*********************************************************/
116
117 VECTOR(head)[0] = 1;
118 for (v = 0; v < no_of_nodes; v++) {
119 VECTOR(next)[v] = v + 2;
120 VECTOR(prev)[v] = v;
121 }
122 VECTOR(next)[no_of_nodes - 1] = 0;
123 /* size is already all zero */
124
125 /***************/
126 /* i:=n; j:=0; */
127 /***************/
128
129 i = no_of_nodes; j = 0;
130
131 /**************/
132 /* do i>=1 -> */
133 /**************/
134
135 while (i >= 1) {
136 long int x, k, len;
137 igraph_vector_int_t *neis;
138
139 /********************************/
140 /* v := delete any from set(j) */
141 /********************************/
142
143 v = VECTOR(head)[j] - 1;
144 x = VECTOR(next)[v];
145 VECTOR(head)[j] = x;
146 if (x != 0) {
147 VECTOR(prev)[x - 1] = 0;
148 }
149
150 /*************************************************/
151 /* alpha(v) := i; alpham1(i) := v; size(v) := -1 */
152 /*************************************************/
153
154 VECTOR(*alpha)[v] = i - 1;
155 if (alpham1) {
156 VECTOR(*alpham1)[i - 1] = v;
157 }
158 VECTOR(size)[v] = -1;
159
160 /********************************************/
161 /* for {v,w} in E such that size(w) >= 0 -> */
162 /********************************************/
163
164 neis = igraph_adjlist_get(&adjlist, v);
165 len = igraph_vector_int_size(neis);
166 for (k = 0; k < len; k++) {
167 long int w = (long int) VECTOR(*neis)[k];
168 long int ws = VECTOR(size)[w];
169 if (ws >= 0) {
170
171 /******************************/
172 /* delete w from set(size(w)) */
173 /******************************/
174
175 long int nw = VECTOR(next)[w];
176 long int pw = VECTOR(prev)[w];
177 if (nw != 0) {
178 VECTOR(prev)[nw - 1] = pw;
179 }
180 if (pw != 0) {
181 VECTOR(next)[pw - 1] = nw;
182 } else {
183 VECTOR(head)[ws] = nw;
184 }
185
186 /******************************/
187 /* size(w) := size(w)+1 */
188 /******************************/
189
190 VECTOR(size)[w] += 1;
191
192 /******************************/
193 /* add w to set(size(w)) */
194 /******************************/
195
196 ws = VECTOR(size)[w];
197 nw = VECTOR(head)[ws];
198 VECTOR(next)[w] = nw;
199 VECTOR(prev)[w] = 0;
200 if (nw != 0) {
201 VECTOR(prev)[nw - 1] = w + 1;
202 }
203 VECTOR(head)[ws] = w + 1;
204
205 }
206 }
207
208 /***********************/
209 /* i := i-1; j := j+1; */
210 /***********************/
211
212 i -= 1;
213 j += 1;
214
215 /*********************************************/
216 /* do j>=0 and set(j)=emptyset -> j:=j-1; od */
217 /*********************************************/
218
219 if (j < no_of_nodes) {
220 while (j >= 0 && VECTOR(head)[j] == 0) {
221 j--;
222 }
223 }
224 }
225
226 igraph_adjlist_destroy(&adjlist);
227 igraph_vector_long_destroy(&prev);
228 igraph_vector_long_destroy(&next);
229 igraph_vector_long_destroy(&head);
230 igraph_vector_long_destroy(&size);
231 IGRAPH_FINALLY_CLEAN(5);
232
233 return 0;
234 }
235
236 /**
237 * \function igraph_is_chordal
238 * Decides whether a graph is chordal
239 *
240 * A graph is chordal if each of its cycles of four or more nodes
241 * has a chord, which is an edge joining two nodes that are not
242 * adjacent in the cycle. An equivalent definition is that any
243 * chordless cycles have at most three nodes.
244 *
245 * If either \p alpha or \p alpha1 is given, then the other is
246 * calculated by taking simply the inverse. If neither are given,
247 * then \ref igraph_maximum_cardinality_search() is called to calculate
248 * them.
249 * \param graph The input graph, it might be directed, but edge
250 * direction is ignored.
251 * \param alpha Either an alpha vector coming from
252 * \ref igraph_maximum_cardinality_search() (on the same graph), or a
253 * null pointer.
254 * \param alpham1 Either an inverse alpha vector coming from \ref
255 * igraph_maximum_cardinality_search() (on the same graph) or a null
256 * pointer.
257 * \param chordal Pointer to a boolean, the result is stored here.
258 * \param fill_in Pointer to an initialized vector, or a null
259 * pointer. If not a null pointer, then the fill-in of the graph is
260 * stored here. The fill-in is the set of edges that are needed to
261 * make the graph chordal. The vector is resized as needed.
262 * \param newgraph Pointer to an uninitialized graph, or a null
263 * pointer. If not a null pointer, then a new triangulated graph is
264 * created here. This essentially means adding the fill-in edges to
265 * the original graph.
266 * \return Error code.
267 *
268 * Time complexity: O(n).
269 *
270 * \sa \ref igraph_maximum_cardinality_search().
271 */
272
igraph_is_chordal(const igraph_t * graph,const igraph_vector_t * alpha,const igraph_vector_t * alpham1,igraph_bool_t * chordal,igraph_vector_t * fill_in,igraph_t * newgraph)273 int igraph_is_chordal(const igraph_t *graph,
274 const igraph_vector_t *alpha,
275 const igraph_vector_t *alpham1,
276 igraph_bool_t *chordal,
277 igraph_vector_t *fill_in,
278 igraph_t *newgraph) {
279
280 long int no_of_nodes = igraph_vcount(graph);
281 const igraph_vector_t *my_alpha = alpha, *my_alpham1 = alpham1;
282 igraph_vector_t v_alpha, v_alpham1;
283 igraph_vector_long_t f, index;
284 long int i;
285 igraph_adjlist_t adjlist;
286 igraph_vector_long_t mark;
287 igraph_bool_t calc_edges = fill_in || newgraph;
288 igraph_vector_t *my_fill_in = fill_in, v_fill_in;
289
290 /*****************/
291 /* local v, w, x */
292 /*****************/
293
294 long int v, w, x;
295
296 if (!chordal && !calc_edges) {
297 /* Nothing to calculate */
298 return 0;
299 }
300
301 if (!alpha && !alpham1) {
302 IGRAPH_VECTOR_INIT_FINALLY(&v_alpha, no_of_nodes);
303 my_alpha = &v_alpha;
304 IGRAPH_VECTOR_INIT_FINALLY(&v_alpham1, no_of_nodes);
305 my_alpham1 = &v_alpham1;
306 IGRAPH_CHECK(igraph_maximum_cardinality_search(graph,
307 (igraph_vector_t*) my_alpha,
308 (igraph_vector_t*) my_alpham1));
309 } else if (alpha && !alpham1) {
310 long int v;
311 IGRAPH_VECTOR_INIT_FINALLY(&v_alpham1, no_of_nodes);
312 my_alpham1 = &v_alpham1;
313 for (v = 0; v < no_of_nodes; v++) {
314 long int i = (long int) VECTOR(*my_alpha)[v];
315 VECTOR(*my_alpham1)[i] = v;
316 }
317 } else if (!alpha && alpham1) {
318 long int i;
319 IGRAPH_VECTOR_INIT_FINALLY(&v_alpha, no_of_nodes);
320 my_alpha = &v_alpha;
321 for (i = 0; i < no_of_nodes; i++) {
322 long int v = (long int) VECTOR(*my_alpham1)[i];
323 VECTOR(*my_alpha)[v] = i;
324 }
325 }
326
327 if (!fill_in && newgraph) {
328 IGRAPH_VECTOR_INIT_FINALLY(&v_fill_in, 0);
329 my_fill_in = &v_fill_in;
330 }
331
332 IGRAPH_CHECK(igraph_vector_long_init(&f, no_of_nodes));
333 IGRAPH_FINALLY(igraph_vector_long_destroy, &f);
334 IGRAPH_CHECK(igraph_vector_long_init(&index, no_of_nodes));
335 IGRAPH_FINALLY(igraph_vector_long_destroy, &index);
336 IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_ALL));
337 IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
338 IGRAPH_CHECK(igraph_vector_long_init(&mark, no_of_nodes));
339 IGRAPH_FINALLY(igraph_vector_long_destroy, &mark);
340 if (my_fill_in) {
341 igraph_vector_clear(my_fill_in);
342 }
343
344 if (chordal) {
345 *chordal = 1;
346 }
347
348 /*********************/
349 /* for i in [1,n] -> */
350 /*********************/
351
352 for (i = 0; i < no_of_nodes; i++) {
353 igraph_vector_int_t *neis;
354 long int j, len;
355
356 /**********************************************/
357 /* w := alpham1(i); f(w) := w; index(w) := i; */
358 /**********************************************/
359
360 w = (long int) VECTOR(*my_alpham1)[i];
361 VECTOR(f)[w] = w;
362 VECTOR(index)[w] = i;
363
364 /******************************************/
365 /* for {v,w} in E such that alpha(v)<i -> */
366 /******************************************/
367
368 neis = igraph_adjlist_get(&adjlist, w);
369 len = igraph_vector_int_size(neis);
370 for (j = 0; j < len; j++) {
371 v = (long int) VECTOR(*neis)[j];
372 VECTOR(mark)[v] = w + 1;
373 }
374
375 for (j = 0; j < len; j++) {
376 v = (long int) VECTOR(*neis)[j];
377 if (VECTOR(*my_alpha)[v] >= i) {
378 continue;
379 }
380
381 /**********/
382 /* x := v */
383 /**********/
384
385 x = v;
386
387 /********************/
388 /* do index(x)<i -> */
389 /********************/
390
391 while (VECTOR(index)[x] < i) {
392
393 /******************/
394 /* index(x) := i; */
395 /******************/
396
397 VECTOR(index)[x] = i;
398
399 /**********************************/
400 /* add {x,w} to E union F(alpha); */
401 /**********************************/
402
403 if (VECTOR(mark)[x] != w + 1) {
404
405 if (chordal) {
406 *chordal = 0;
407 }
408
409 if (my_fill_in) {
410 IGRAPH_CHECK(igraph_vector_push_back(my_fill_in, x));
411 IGRAPH_CHECK(igraph_vector_push_back(my_fill_in, w));
412 }
413
414 if (!calc_edges) {
415 /* make sure that we exit from all loops */
416 i = no_of_nodes;
417 j = len;
418 break;
419 }
420 }
421
422 /*************/
423 /* x := f(x) */
424 /*************/
425
426 x = VECTOR(f)[x];
427
428 } /* while (VECTOR(index)[x] < i) */
429
430 /*****************************/
431 /* if (f(x)=x -> f(x):=w; fi */
432 /*****************************/
433
434 if (VECTOR(f)[x] == x) {
435 VECTOR(f)[x] = w;
436 }
437 }
438 }
439
440 igraph_vector_long_destroy(&mark);
441 igraph_adjlist_destroy(&adjlist);
442 igraph_vector_long_destroy(&index);
443 igraph_vector_long_destroy(&f);
444 IGRAPH_FINALLY_CLEAN(4);
445
446 if (newgraph) {
447 IGRAPH_CHECK(igraph_copy(newgraph, graph));
448 IGRAPH_FINALLY(igraph_destroy, newgraph);
449 IGRAPH_CHECK(igraph_add_edges(newgraph, my_fill_in, 0));
450 IGRAPH_FINALLY_CLEAN(1);
451 }
452
453 if (!fill_in && newgraph) {
454 igraph_vector_destroy(&v_fill_in);
455 IGRAPH_FINALLY_CLEAN(1);
456 }
457
458 if (!alpha && !alpham1) {
459 igraph_vector_destroy(&v_alpham1);
460 igraph_vector_destroy(&v_alpha);
461 IGRAPH_FINALLY_CLEAN(2);
462 } else if (alpha && !alpham1) {
463 igraph_vector_destroy(&v_alpham1);
464 IGRAPH_FINALLY_CLEAN(1);
465 } else if (!alpha && alpham1) {
466 igraph_vector_destroy(&v_alpha);
467 IGRAPH_FINALLY_CLEAN(1);
468 }
469
470 return 0;
471 }
472