1 /* -*- mode: C -*- */
2 /*
3 IGraph library.
4 Copyright (C) 2003-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_components.h"
25
26 #include "igraph_adjlist.h"
27 #include "igraph_dqueue.h"
28 #include "igraph_interface.h"
29 #include "igraph_memory.h"
30 #include "igraph_operators.h"
31 #include "igraph_progress.h"
32 #include "igraph_stack.h"
33 #include "igraph_structural.h"
34 #include "igraph_vector.h"
35
36 #include "core/interruption.h"
37 #include "operators/subgraph.h"
38
39 #include <limits.h>
40
41 static int igraph_i_clusters_weak(const igraph_t *graph, igraph_vector_t *membership,
42 igraph_vector_t *csize, igraph_integer_t *no);
43
44 static int igraph_i_clusters_strong(const igraph_t *graph, igraph_vector_t *membership,
45 igraph_vector_t *csize, igraph_integer_t *no);
46
47 /**
48 * \ingroup structural
49 * \function igraph_clusters
50 * \brief Calculates the (weakly or strongly) connected components in a graph.
51 *
52 * \param graph The graph object to analyze.
53 * \param membership First half of the result will be stored here. For
54 * every vertex the id of its component is given. The vector
55 * has to be preinitialized and will be resized. Alternatively
56 * this argument can be \c NULL, in which case it is ignored.
57 * \param csize The second half of the result. For every component it
58 * gives its size, the order is defined by the component ids.
59 * The vector has to be preinitialized and will be resized.
60 * Alternatively this argument can be \c NULL, in which
61 * case it is ignored.
62 * \param no Pointer to an integer, if not \c NULL then the number of
63 * clusters will be stored here.
64 * \param mode For directed graph this specifies whether to calculate
65 * weakly or strongly connected components. Possible values:
66 * \c IGRAPH_WEAK,
67 * \c IGRAPH_STRONG. This argument is
68 * ignored for undirected graphs.
69 * \return Error code:
70 * \c IGRAPH_EINVAL: invalid mode argument.
71 *
72 * Time complexity: O(|V|+|E|),
73 * |V| and
74 * |E| are the number of vertices and
75 * edges in the graph.
76 */
77
igraph_clusters(const igraph_t * graph,igraph_vector_t * membership,igraph_vector_t * csize,igraph_integer_t * no,igraph_connectedness_t mode)78 int igraph_clusters(const igraph_t *graph, igraph_vector_t *membership,
79 igraph_vector_t *csize, igraph_integer_t *no,
80 igraph_connectedness_t mode) {
81 if (mode == IGRAPH_WEAK || !igraph_is_directed(graph)) {
82 return igraph_i_clusters_weak(graph, membership, csize, no);
83 } else if (mode == IGRAPH_STRONG) {
84 return igraph_i_clusters_strong(graph, membership, csize, no);
85 }
86
87 IGRAPH_ERROR("Cannot calculate clusters", IGRAPH_EINVAL);
88 }
89
igraph_i_clusters_weak(const igraph_t * graph,igraph_vector_t * membership,igraph_vector_t * csize,igraph_integer_t * no)90 static int igraph_i_clusters_weak(const igraph_t *graph, igraph_vector_t *membership,
91 igraph_vector_t *csize, igraph_integer_t *no) {
92
93 long int no_of_nodes = igraph_vcount(graph);
94 char *already_added;
95 long int first_node, act_cluster_size = 0, no_of_clusters = 1;
96
97 igraph_dqueue_t q = IGRAPH_DQUEUE_NULL;
98
99 long int i;
100 igraph_vector_t neis = IGRAPH_VECTOR_NULL;
101
102 already_added = IGRAPH_CALLOC(no_of_nodes, char);
103 if (already_added == 0) {
104 IGRAPH_ERROR("Cannot calculate clusters", IGRAPH_ENOMEM);
105 }
106 IGRAPH_FINALLY(igraph_free, already_added);
107
108 IGRAPH_DQUEUE_INIT_FINALLY(&q, no_of_nodes > 100000 ? 10000 : no_of_nodes / 10);
109 IGRAPH_VECTOR_INIT_FINALLY(&neis, 0);
110
111 /* Memory for result, csize is dynamically allocated */
112 if (membership) {
113 IGRAPH_CHECK(igraph_vector_resize(membership, no_of_nodes));
114 }
115 if (csize) {
116 igraph_vector_clear(csize);
117 }
118
119 /* The algorithm */
120
121 for (first_node = 0; first_node < no_of_nodes; ++first_node) {
122 if (already_added[first_node] == 1) {
123 continue;
124 }
125 IGRAPH_ALLOW_INTERRUPTION();
126
127 already_added[first_node] = 1;
128 act_cluster_size = 1;
129 if (membership) {
130 VECTOR(*membership)[first_node] = no_of_clusters - 1;
131 }
132 IGRAPH_CHECK(igraph_dqueue_push(&q, first_node));
133
134 while ( !igraph_dqueue_empty(&q) ) {
135 long int act_node = (long int) igraph_dqueue_pop(&q);
136 IGRAPH_CHECK(igraph_neighbors(graph, &neis,
137 (igraph_integer_t) act_node, IGRAPH_ALL));
138 for (i = 0; i < igraph_vector_size(&neis); i++) {
139 long int neighbor = (long int) VECTOR(neis)[i];
140 if (already_added[neighbor] == 1) {
141 continue;
142 }
143 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
144 already_added[neighbor] = 1;
145 act_cluster_size++;
146 if (membership) {
147 VECTOR(*membership)[neighbor] = no_of_clusters - 1;
148 }
149 }
150 }
151 no_of_clusters++;
152 if (csize) {
153 IGRAPH_CHECK(igraph_vector_push_back(csize, act_cluster_size));
154 }
155 }
156
157 /* Cleaning up */
158
159 if (no) {
160 *no = (igraph_integer_t) no_of_clusters - 1;
161 }
162
163 IGRAPH_FREE(already_added);
164 igraph_dqueue_destroy(&q);
165 igraph_vector_destroy(&neis);
166 IGRAPH_FINALLY_CLEAN(3);
167
168 return 0;
169 }
170
igraph_i_clusters_strong(const igraph_t * graph,igraph_vector_t * membership,igraph_vector_t * csize,igraph_integer_t * no)171 static int igraph_i_clusters_strong(const igraph_t *graph, igraph_vector_t *membership,
172 igraph_vector_t *csize, igraph_integer_t *no) {
173
174 long int no_of_nodes = igraph_vcount(graph);
175 igraph_vector_t next_nei = IGRAPH_VECTOR_NULL;
176
177 long int i, n, num_seen;
178 igraph_dqueue_t q = IGRAPH_DQUEUE_NULL;
179
180 long int no_of_clusters = 1;
181 long int act_cluster_size;
182
183 igraph_vector_t out = IGRAPH_VECTOR_NULL;
184 const igraph_vector_int_t* tmp;
185
186 igraph_adjlist_t adjlist;
187
188 /* The result */
189
190 IGRAPH_VECTOR_INIT_FINALLY(&next_nei, no_of_nodes);
191 IGRAPH_VECTOR_INIT_FINALLY(&out, 0);
192 IGRAPH_DQUEUE_INIT_FINALLY(&q, 100);
193
194 if (membership) {
195 IGRAPH_CHECK(igraph_vector_resize(membership, no_of_nodes));
196 }
197 IGRAPH_CHECK(igraph_vector_reserve(&out, no_of_nodes));
198
199 igraph_vector_null(&out);
200 if (csize) {
201 igraph_vector_clear(csize);
202 }
203
204 IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_OUT, IGRAPH_LOOPS_ONCE, IGRAPH_MULTIPLE));
205 IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
206
207 num_seen = 0;
208 for (i = 0; i < no_of_nodes; i++) {
209 IGRAPH_ALLOW_INTERRUPTION();
210
211 tmp = igraph_adjlist_get(&adjlist, i);
212 if (VECTOR(next_nei)[i] > igraph_vector_int_size(tmp)) {
213 continue;
214 }
215
216 IGRAPH_CHECK(igraph_dqueue_push(&q, i));
217 while (!igraph_dqueue_empty(&q)) {
218 long int act_node = (long int) igraph_dqueue_back(&q);
219 tmp = igraph_adjlist_get(&adjlist, act_node);
220 if (VECTOR(next_nei)[act_node] == 0) {
221 /* this is the first time we've met this vertex */
222 VECTOR(next_nei)[act_node]++;
223 } else if (VECTOR(next_nei)[act_node] <= igraph_vector_int_size(tmp)) {
224 /* we've already met this vertex but it has more children */
225 long int neighbor = (long int) VECTOR(*tmp)[(long int)
226 VECTOR(next_nei)[act_node] - 1];
227 if (VECTOR(next_nei)[neighbor] == 0) {
228 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
229 }
230 VECTOR(next_nei)[act_node]++;
231 } else {
232 /* we've met this vertex and it has no more children */
233 IGRAPH_CHECK(igraph_vector_push_back(&out, act_node));
234 igraph_dqueue_pop_back(&q);
235 num_seen++;
236
237 if (num_seen % 10000 == 0) {
238 /* time to report progress and allow the user to interrupt */
239 IGRAPH_PROGRESS("Strongly connected components: ",
240 num_seen * 50.0 / no_of_nodes, NULL);
241 IGRAPH_ALLOW_INTERRUPTION();
242 }
243 }
244 } /* while q */
245 } /* for */
246
247 IGRAPH_PROGRESS("Strongly connected components: ", 50.0, NULL);
248
249 igraph_adjlist_destroy(&adjlist);
250 IGRAPH_FINALLY_CLEAN(1);
251
252 IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_IN, IGRAPH_LOOPS_ONCE, IGRAPH_MULTIPLE));
253 IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
254
255 /* OK, we've the 'out' values for the nodes, let's use them in
256 decreasing order with the help of a heap */
257
258 igraph_vector_null(&next_nei); /* mark already added vertices */
259 num_seen = 0;
260
261 while (!igraph_vector_empty(&out)) {
262 long int grandfather = (long int) igraph_vector_pop_back(&out);
263
264 if (VECTOR(next_nei)[grandfather] != 0) {
265 continue;
266 }
267 VECTOR(next_nei)[grandfather] = 1;
268 act_cluster_size = 1;
269 if (membership) {
270 VECTOR(*membership)[grandfather] = no_of_clusters - 1;
271 }
272 IGRAPH_CHECK(igraph_dqueue_push(&q, grandfather));
273
274 num_seen++;
275 if (num_seen % 10000 == 0) {
276 /* time to report progress and allow the user to interrupt */
277 IGRAPH_PROGRESS("Strongly connected components: ",
278 50.0 + num_seen * 50.0 / no_of_nodes, NULL);
279 IGRAPH_ALLOW_INTERRUPTION();
280 }
281
282 while (!igraph_dqueue_empty(&q)) {
283 long int act_node = (long int) igraph_dqueue_pop_back(&q);
284 tmp = igraph_adjlist_get(&adjlist, act_node);
285 n = igraph_vector_int_size(tmp);
286 for (i = 0; i < n; i++) {
287 long int neighbor = (long int) VECTOR(*tmp)[i];
288 if (VECTOR(next_nei)[neighbor] != 0) {
289 continue;
290 }
291 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
292 VECTOR(next_nei)[neighbor] = 1;
293 act_cluster_size++;
294 if (membership) {
295 VECTOR(*membership)[neighbor] = no_of_clusters - 1;
296 }
297
298 num_seen++;
299 if (num_seen % 10000 == 0) {
300 /* time to report progress and allow the user to interrupt */
301 IGRAPH_PROGRESS("Strongly connected components: ",
302 50.0 + num_seen * 50.0 / no_of_nodes, NULL);
303 IGRAPH_ALLOW_INTERRUPTION();
304 }
305 }
306 }
307
308 no_of_clusters++;
309 if (csize) {
310 IGRAPH_CHECK(igraph_vector_push_back(csize, act_cluster_size));
311 }
312 }
313
314 IGRAPH_PROGRESS("Strongly connected components: ", 100.0, NULL);
315
316 if (no) {
317 *no = (igraph_integer_t) no_of_clusters - 1;
318 }
319
320 /* Clean up, return */
321
322 igraph_adjlist_destroy(&adjlist);
323 igraph_vector_destroy(&out);
324 igraph_dqueue_destroy(&q);
325 igraph_vector_destroy(&next_nei);
326 IGRAPH_FINALLY_CLEAN(4);
327
328 return 0;
329 }
330
331 int igraph_is_connected_weak(const igraph_t *graph, igraph_bool_t *res);
332
333 /**
334 * \ingroup structural
335 * \function igraph_is_connected
336 * \brief Decides whether the graph is (weakly or strongly) connected.
337 *
338 * A graph with zero vertices (i.e. the null graph) is \em not connected by
339 * definition. This behaviour changed in igraph 0.9; earlier versions assumed
340 * that the null graph is connected. See the following issue on Github for the
341 * argument that led us to change the definition:
342 * https://github.com/igraph/igraph/issues/1538
343 *
344 * \param graph The graph object to analyze.
345 * \param res Pointer to a logical variable, the result will be stored
346 * here.
347 * \param mode For a directed graph this specifies whether to calculate
348 * weak or strong connectedness. Possible values:
349 * \c IGRAPH_WEAK,
350 * \c IGRAPH_STRONG. This argument is
351 * ignored for undirected graphs.
352 * \return Error code:
353 * \c IGRAPH_EINVAL: invalid mode argument.
354 *
355 * Time complexity: O(|V|+|E|), the
356 * number of vertices
357 * plus the number of edges in the graph.
358 */
359
igraph_is_connected(const igraph_t * graph,igraph_bool_t * res,igraph_connectedness_t mode)360 int igraph_is_connected(const igraph_t *graph, igraph_bool_t *res,
361 igraph_connectedness_t mode) {
362 if (igraph_vcount(graph) == 0) {
363 /* Changed in igraph 0.9; see https://github.com/igraph/igraph/issues/1538
364 * for the reasoning behind the change */
365 *res = 0;
366 return IGRAPH_SUCCESS;
367 }
368
369 if (mode == IGRAPH_WEAK || !igraph_is_directed(graph)) {
370 return igraph_is_connected_weak(graph, res);
371 } else if (mode == IGRAPH_STRONG) {
372 int retval;
373 igraph_integer_t no;
374 retval = igraph_i_clusters_strong(graph, 0, 0, &no);
375 *res = (no == 1);
376 return retval;
377 }
378
379 IGRAPH_ERROR("mode argument", IGRAPH_EINVAL);
380 }
381
382 /**
383 * \ingroup structural
384 * \function igraph_is_connected_weak
385 * \brief Query whether the graph is weakly connected.
386 *
387 * A graph with zero vertices (i.e. the null graph) is weakly connected by
388 * definition. A directed graph is weakly connected if its undirected version
389 * is connected. In the case of undirected graphs, weakly connected and
390 * connected are equivalent.
391 *
392 * \param graph The graph object to analyze.
393 * \param res Pointer to a logical variable; the result will be stored here.
394 * \return Error code:
395 * \c IGRAPH_ENOMEM: unable to allocate requested memory.
396 *
397 * Time complexity: O(|V|+|E|), the number of vertices plus the number of
398 * edges in the graph.
399 */
400
igraph_is_connected_weak(const igraph_t * graph,igraph_bool_t * res)401 int igraph_is_connected_weak(const igraph_t *graph, igraph_bool_t *res) {
402
403 long int no_of_nodes = igraph_vcount(graph);
404 char *already_added;
405 igraph_vector_t neis = IGRAPH_VECTOR_NULL;
406 igraph_dqueue_t q = IGRAPH_DQUEUE_NULL;
407
408 long int i, j;
409
410 if (no_of_nodes == 0) {
411 *res = 1;
412 return IGRAPH_SUCCESS;
413 }
414
415 already_added = IGRAPH_CALLOC(no_of_nodes, char);
416 if (already_added == 0) {
417 IGRAPH_ERROR("is connected (weak) failed", IGRAPH_ENOMEM);
418 }
419 IGRAPH_FINALLY(igraph_free, already_added);
420
421 IGRAPH_DQUEUE_INIT_FINALLY(&q, 10);
422 IGRAPH_VECTOR_INIT_FINALLY(&neis, 0);
423
424 /* Try to find at least two clusters */
425 already_added[0] = 1;
426 IGRAPH_CHECK(igraph_dqueue_push(&q, 0));
427
428 j = 1;
429 while ( !igraph_dqueue_empty(&q)) {
430 long int actnode = (long int) igraph_dqueue_pop(&q);
431 IGRAPH_ALLOW_INTERRUPTION();
432 IGRAPH_CHECK(igraph_neighbors(graph, &neis, (igraph_integer_t) actnode,
433 IGRAPH_ALL));
434 for (i = 0; i < igraph_vector_size(&neis); i++) {
435 long int neighbor = (long int) VECTOR(neis)[i];
436 if (already_added[neighbor] != 0) {
437 continue;
438 }
439 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
440 j++;
441 already_added[neighbor]++;
442 }
443 }
444
445 /* Connected? */
446 *res = (j == no_of_nodes);
447
448 IGRAPH_FREE(already_added);
449 igraph_dqueue_destroy(&q);
450 igraph_vector_destroy(&neis);
451 IGRAPH_FINALLY_CLEAN(3);
452
453 return 0;
454 }
455
456 /**
457 * \function igraph_decompose_destroy
458 * \brief Free the memory allocated by \ref igraph_decompose().
459 *
460 * \param complist The list of graph components, as returned by
461 * \ref igraph_decompose().
462 *
463 * Time complexity: O(c), c is the number of components.
464 */
465
igraph_decompose_destroy(igraph_vector_ptr_t * complist)466 void igraph_decompose_destroy(igraph_vector_ptr_t *complist) {
467 long int i;
468 for (i = 0; i < igraph_vector_ptr_size(complist); i++) {
469 if (VECTOR(*complist)[i] != 0) {
470 igraph_destroy(VECTOR(*complist)[i]);
471 igraph_free(VECTOR(*complist)[i]);
472 }
473 }
474 }
475
476 static int igraph_i_decompose_weak(const igraph_t *graph,
477 igraph_vector_ptr_t *components,
478 long int maxcompno, long int minelements);
479
480 static int igraph_i_decompose_strong(const igraph_t *graph,
481 igraph_vector_ptr_t *components,
482 long int maxcompno, long int minelements);
483
484 /**
485 * \function igraph_decompose
486 * \brief Decompose a graph into connected components.
487 *
488 * Create separate graph for each component of a graph. Note that the
489 * vertex ids in the new graphs will be different than in the original
490 * graph. (Except if there is only one component in the original graph.)
491 *
492 * \param graph The original graph.
493 * \param components This pointer vector will contain pointers to the
494 * subcomponent graphs. It should be initialized before calling this
495 * function and will be resized to hold the graphs. Don't forget to
496 * call \ref igraph_destroy() and \ref igraph_free() on the elements of
497 * this pointer vector to free unneeded memory. Alternatively, you can
498 * simply call \ref igraph_decompose_destroy() that does this for you.
499 * \param mode Either \c IGRAPH_WEAK or \c IGRAPH_STRONG for weakly
500 * and strongly connected components respectively.
501 * \param maxcompno The maximum number of components to return. The
502 * first \p maxcompno components will be returned (which hold at
503 * least \p minelements vertices, see the next parameter), the
504 * others will be ignored. Supply -1 here if you don't want to limit
505 * the number of components.
506 * \param minelements The minimum number of vertices a component
507 * should contain in order to place it in the \p components
508 * vector. Eg. supply 2 here to ignore isolated vertices.
509 * \return Error code, \c IGRAPH_ENOMEM if there is not enough memory
510 * to perform the operation.
511 *
512 * Added in version 0.2.</para><para>
513 *
514 * Time complexity: O(|V|+|E|), the number of vertices plus the number
515 * of edges.
516 *
517 * \example examples/simple/igraph_decompose.c
518 */
519
igraph_decompose(const igraph_t * graph,igraph_vector_ptr_t * components,igraph_connectedness_t mode,long int maxcompno,long int minelements)520 int igraph_decompose(const igraph_t *graph, igraph_vector_ptr_t *components,
521 igraph_connectedness_t mode,
522 long int maxcompno, long int minelements) {
523 if (mode == IGRAPH_WEAK || !igraph_is_directed(graph)) {
524 return igraph_i_decompose_weak(graph, components, maxcompno, minelements);
525 } else if (mode == IGRAPH_STRONG) {
526 return igraph_i_decompose_strong(graph, components, maxcompno, minelements);
527 }
528
529 IGRAPH_ERROR("Cannot decompose graph", IGRAPH_EINVAL);
530 }
531
igraph_i_decompose_weak(const igraph_t * graph,igraph_vector_ptr_t * components,long int maxcompno,long int minelements)532 static int igraph_i_decompose_weak(const igraph_t *graph,
533 igraph_vector_ptr_t *components,
534 long int maxcompno, long int minelements) {
535
536 long int actstart;
537 long int no_of_nodes = igraph_vcount(graph);
538 long int resco = 0; /* number of graphs created so far */
539 char *already_added;
540 igraph_dqueue_t q;
541 igraph_vector_t verts;
542 igraph_vector_t neis;
543 igraph_vector_t vids_old2new;
544 long int i;
545 igraph_t *newg;
546
547
548 if (maxcompno < 0) {
549 maxcompno = LONG_MAX;
550 }
551
552 igraph_vector_ptr_clear(components);
553 IGRAPH_FINALLY(igraph_decompose_destroy, components);
554
555 /* already_added keeps track of what nodes made it into a graph already */
556 already_added = IGRAPH_CALLOC(no_of_nodes, char);
557 if (already_added == 0) {
558 IGRAPH_ERROR("Cannot decompose graph", IGRAPH_ENOMEM);
559 }
560 IGRAPH_FINALLY(igraph_free, already_added);
561
562 IGRAPH_CHECK(igraph_dqueue_init(&q, 100));
563 IGRAPH_FINALLY(igraph_dqueue_destroy, &q);
564 IGRAPH_VECTOR_INIT_FINALLY(&verts, 0);
565 IGRAPH_VECTOR_INIT_FINALLY(&neis, 0);
566 IGRAPH_VECTOR_INIT_FINALLY(&vids_old2new, no_of_nodes);
567
568 /* vids_old2new would have been created internally in igraph_induced_subgraph(),
569 but it is slow if the graph is large and consists of many small components,
570 so we create it once here and then re-use it */
571
572 /* add a node and its neighbors at once, recursively
573 then switch to next node that has not been added already */
574 for (actstart = 0; resco < maxcompno && actstart < no_of_nodes; actstart++) {
575
576 if (already_added[actstart]) {
577 continue;
578 }
579 IGRAPH_ALLOW_INTERRUPTION();
580
581 igraph_vector_clear(&verts);
582
583 /* add the node itself */
584 already_added[actstart] = 1;
585 IGRAPH_CHECK(igraph_vector_push_back(&verts, actstart));
586 IGRAPH_CHECK(igraph_dqueue_push(&q, actstart));
587
588 /* add the neighbors, recursively */
589 while (!igraph_dqueue_empty(&q) ) {
590 /* pop from the queue of this component */
591 long int actvert = (long int) igraph_dqueue_pop(&q);
592 IGRAPH_CHECK(igraph_neighbors(graph, &neis, (igraph_integer_t) actvert,
593 IGRAPH_ALL));
594 /* iterate over the neighbors */
595 for (i = 0; i < igraph_vector_size(&neis); i++) {
596 long int neighbor = (long int) VECTOR(neis)[i];
597 if (already_added[neighbor] == 1) {
598 continue;
599 }
600 /* add neighbor */
601 already_added[neighbor] = 1;
602
603 /* recursion: append neighbor to the queues */
604 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
605 IGRAPH_CHECK(igraph_vector_push_back(&verts, neighbor));
606 }
607 }
608
609 /* ok, we have a component */
610 if (igraph_vector_size(&verts) < minelements) {
611 continue;
612 }
613
614 newg = IGRAPH_CALLOC(1, igraph_t);
615 if (newg == 0) {
616 IGRAPH_ERROR("Cannot decompose graph", IGRAPH_ENOMEM);
617 }
618 IGRAPH_CHECK(igraph_vector_ptr_push_back(components, newg));
619 IGRAPH_CHECK(igraph_i_induced_subgraph_map(
620 graph, newg, igraph_vss_vector(&verts),
621 IGRAPH_SUBGRAPH_AUTO, &vids_old2new,
622 /* invmap = */ 0, /* map_is_prepared = */ 1
623 ));
624 resco++;
625
626 /* vids_old2new does not have to be cleaned up here; since we are doing
627 * weak decomposition, each vertex will appear in only one of the
628 * connected components so we won't ever touch an item in vids_old2new
629 * if it was already set to a non-zero value in a previous component */
630
631 } /* for actstart++ */
632
633 igraph_vector_destroy(&vids_old2new);
634 igraph_vector_destroy(&neis);
635 igraph_vector_destroy(&verts);
636 igraph_dqueue_destroy(&q);
637 IGRAPH_FREE(already_added);
638 IGRAPH_FINALLY_CLEAN(6); /* + components */
639
640 return 0;
641 }
642
igraph_i_decompose_strong(const igraph_t * graph,igraph_vector_ptr_t * components,long int maxcompno,long int minelements)643 static int igraph_i_decompose_strong(const igraph_t *graph,
644 igraph_vector_ptr_t *components,
645 long int maxcompno, long int minelements) {
646
647
648 long int no_of_nodes = igraph_vcount(graph);
649
650 /* this is a heap used twice for checking what nodes have
651 * been counted already */
652 igraph_vector_t next_nei = IGRAPH_VECTOR_NULL;
653
654 long int i, n, num_seen;
655 igraph_dqueue_t q = IGRAPH_DQUEUE_NULL;
656
657 long int no_of_clusters = 0;
658 long int act_cluster_size;
659
660 igraph_vector_t out = IGRAPH_VECTOR_NULL;
661 const igraph_vector_int_t* tmp;
662
663 igraph_adjlist_t adjlist;
664 igraph_vector_t verts;
665 igraph_vector_t vids_old2new;
666 igraph_t *newg;
667
668 if (maxcompno < 0) {
669 maxcompno = LONG_MAX;
670 }
671
672 igraph_vector_ptr_clear(components);
673 IGRAPH_FINALLY(igraph_decompose_destroy, components);
674
675 /* The result */
676
677 IGRAPH_VECTOR_INIT_FINALLY(&vids_old2new, no_of_nodes);
678 IGRAPH_VECTOR_INIT_FINALLY(&verts, 0);
679 IGRAPH_VECTOR_INIT_FINALLY(&next_nei, no_of_nodes);
680 IGRAPH_VECTOR_INIT_FINALLY(&out, 0);
681 IGRAPH_DQUEUE_INIT_FINALLY(&q, 100);
682
683 IGRAPH_CHECK(igraph_vector_reserve(&out, no_of_nodes));
684
685 igraph_vector_null(&out);
686
687 IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_OUT, IGRAPH_LOOPS_ONCE, IGRAPH_MULTIPLE));
688 IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
689
690 /* vids_old2new would have been created internally in igraph_induced_subgraph(),
691 but it is slow if the graph is large and consists of many small components,
692 so we create it once here and then re-use it */
693
694 /* number of components seen */
695 num_seen = 0;
696 /* populate the 'out' vector by browsing a node and following up
697 all its neighbors recursively, then switching to the next
698 unassigned node */
699 for (i = 0; i < no_of_nodes; i++) {
700 IGRAPH_ALLOW_INTERRUPTION();
701
702 /* get all the 'out' neighbors of this node
703 * NOTE: next_nei is initialized [0, 0, ...] */
704 tmp = igraph_adjlist_get(&adjlist, i);
705 if (VECTOR(next_nei)[i] > igraph_vector_int_size(tmp)) {
706 continue;
707 }
708
709 /* add this node to the queue for this component */
710 IGRAPH_CHECK(igraph_dqueue_push(&q, i));
711
712 /* consume the tree from this node ("root") recursively
713 * until there is no more */
714 while (!igraph_dqueue_empty(&q)) {
715 /* this looks up but does NOT consume the queue */
716 long int act_node = (long int) igraph_dqueue_back(&q);
717
718 /* get all neighbors of this node */
719 tmp = igraph_adjlist_get(&adjlist, act_node);
720 if (VECTOR(next_nei)[act_node] == 0) {
721 /* this is the first time we've met this vertex,
722 * because next_nei is initialized [0, 0, ...] */
723 VECTOR(next_nei)[act_node]++;
724 /* back to the queue, same vertex is up again */
725
726 } else if (VECTOR(next_nei)[act_node] <= igraph_vector_int_size(tmp)) {
727 /* we've already met this vertex but it has more children */
728 long int neighbor = (long int) VECTOR(*tmp)[(long int)
729 VECTOR(next_nei)[act_node] - 1];
730 if (VECTOR(next_nei)[neighbor] == 0) {
731 /* add the root of the other children to the queue */
732 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
733 }
734 VECTOR(next_nei)[act_node]++;
735 } else {
736 /* we've met this vertex and it has no more children */
737 IGRAPH_CHECK(igraph_vector_push_back(&out, act_node));
738 /* this consumes the queue, since there's nowhere to go */
739 igraph_dqueue_pop_back(&q);
740 num_seen++;
741
742 if (num_seen % 10000 == 0) {
743 /* time to report progress and allow the user to interrupt */
744 IGRAPH_PROGRESS("Strongly connected components: ",
745 num_seen * 50.0 / no_of_nodes, NULL);
746 IGRAPH_ALLOW_INTERRUPTION();
747 }
748 }
749 } /* while q */
750 } /* for */
751
752 IGRAPH_PROGRESS("Strongly connected components: ", 50.0, NULL);
753
754 igraph_adjlist_destroy(&adjlist);
755 IGRAPH_FINALLY_CLEAN(1);
756
757 IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_IN, IGRAPH_LOOPS_ONCE, IGRAPH_MULTIPLE));
758 IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
759
760 /* OK, we've the 'out' values for the nodes, let's use them in
761 * decreasing order with the help of the next_nei heap */
762
763 igraph_vector_null(&next_nei); /* mark already added vertices */
764
765 /* number of components built */
766 num_seen = 0;
767 while (!igraph_vector_empty(&out) && no_of_clusters < maxcompno) {
768 /* consume the vector from the last element */
769 long int grandfather = (long int) igraph_vector_pop_back(&out);
770
771 /* been here, done that
772 * NOTE: next_nei is initialized as [0, 0, ...] */
773 if (VECTOR(next_nei)[grandfather] != 0) {
774 continue;
775 }
776
777 /* collect all the members of this component */
778 igraph_vector_clear(&verts);
779
780 /* this node is gone for any future components */
781 VECTOR(next_nei)[grandfather] = 1;
782 act_cluster_size = 1;
783
784 /* add to component */
785 IGRAPH_CHECK(igraph_vector_push_back(&verts, grandfather));
786 IGRAPH_CHECK(igraph_dqueue_push(&q, grandfather));
787
788 num_seen++;
789 if (num_seen % 10000 == 0) {
790 /* time to report progress and allow the user to interrupt */
791 IGRAPH_PROGRESS("Strongly connected components: ",
792 50.0 + num_seen * 50.0 / no_of_nodes, NULL);
793 IGRAPH_ALLOW_INTERRUPTION();
794 }
795
796 while (!igraph_dqueue_empty(&q)) {
797 /* consume the queue from this node */
798 long int act_node = (long int) igraph_dqueue_pop_back(&q);
799 tmp = igraph_adjlist_get(&adjlist, act_node);
800 n = igraph_vector_int_size(tmp);
801 for (i = 0; i < n; i++) {
802 long int neighbor = (long int) VECTOR(*tmp)[i];
803 if (VECTOR(next_nei)[neighbor] != 0) {
804 continue;
805 }
806 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
807 VECTOR(next_nei)[neighbor] = 1;
808 act_cluster_size++;
809
810 /* add to component */
811 IGRAPH_CHECK(igraph_vector_push_back(&verts, neighbor));
812
813 num_seen++;
814 if (num_seen % 10000 == 0) {
815 /* time to report progress and allow the user to interrupt */
816 IGRAPH_PROGRESS("Strongly connected components: ",
817 50.0 + num_seen * 50.0 / no_of_nodes, NULL);
818 IGRAPH_ALLOW_INTERRUPTION();
819 }
820 }
821 }
822
823 /* ok, we have a component */
824 if (igraph_vector_size(&verts) < minelements) {
825 continue;
826 }
827
828 newg = IGRAPH_CALLOC(1, igraph_t);
829 if (newg == 0) {
830 IGRAPH_ERROR("Cannot decompose graph", IGRAPH_ENOMEM);
831 }
832 IGRAPH_CHECK(igraph_vector_ptr_push_back(components, newg));
833 IGRAPH_CHECK(igraph_i_induced_subgraph_map(
834 graph, newg, igraph_vss_vector(&verts),
835 IGRAPH_SUBGRAPH_AUTO, &vids_old2new,
836 /* invmap = */ 0, /* map_is_prepared = */ 1
837 ));
838
839 /* vids_old2new has to be cleaned up here because a vertex may appear
840 * in multiple strongly connected components. Simply calling
841 * igraph_vector_fill() would be an O(n) operation where n is the number
842 * of vertices in the large graph so we cannot do that; we have to
843 * iterate over 'verts' instead */
844 n = igraph_vector_size(&verts);
845 for (i = 0; i < n; i++) {
846 VECTOR(vids_old2new)[(igraph_integer_t) VECTOR(verts)[i]] = 0;
847 }
848
849 no_of_clusters++;
850 }
851
852 IGRAPH_PROGRESS("Strongly connected components: ", 100.0, NULL);
853
854 /* Clean up, return */
855
856 igraph_vector_destroy(&vids_old2new);
857 igraph_vector_destroy(&verts);
858 igraph_adjlist_destroy(&adjlist);
859 igraph_vector_destroy(&out);
860 igraph_dqueue_destroy(&q);
861 igraph_vector_destroy(&next_nei);
862 IGRAPH_FINALLY_CLEAN(7); /* + components */
863
864 return 0;
865
866 }
867
868 /**
869 * \function igraph_articulation_points
870 * Find the articulation points in a graph.
871 *
872 * A vertex is an articulation point if its removal increases
873 * the number of connected components in the graph.
874 * \param graph The input graph.
875 * \param res Pointer to an initialized vector, the
876 * articulation points will be stored here.
877 * \return Error code.
878 *
879 * Time complexity: O(|V|+|E|), linear in the number of vertices and edges.
880 *
881 * \sa \ref igraph_biconnected_components(), \ref igraph_clusters(), \ref igraph_bridges()
882 */
883
igraph_articulation_points(const igraph_t * graph,igraph_vector_t * res)884 int igraph_articulation_points(const igraph_t *graph,
885 igraph_vector_t *res) {
886
887 igraph_integer_t no;
888 return igraph_biconnected_components(graph, &no, 0, 0, 0, res);
889 }
890
891 void igraph_i_free_vectorlist(igraph_vector_ptr_t *list);
892
igraph_i_free_vectorlist(igraph_vector_ptr_t * list)893 void igraph_i_free_vectorlist(igraph_vector_ptr_t *list) {
894 long int i, n = igraph_vector_ptr_size(list);
895 for (i = 0; i < n; i++) {
896 igraph_vector_t *v = VECTOR(*list)[i];
897 if (v) {
898 igraph_vector_destroy(v);
899 IGRAPH_FREE(v);
900 }
901 }
902 igraph_vector_ptr_destroy(list);
903 }
904
905 /**
906 * \function igraph_biconnected_components
907 * Calculate biconnected components
908 *
909 * A graph is biconnected if the removal of any single vertex (and
910 * its incident edges) does not disconnect it.
911 *
912 * </para><para>
913 * A biconnected component of a graph is a maximal biconnected
914 * subgraph of it. The biconnected components of a graph can be given
915 * by the partition of its edges: every edge is a member of exactly
916 * one biconnected component. Note that this is not true for
917 * vertices: the same vertex can be part of many biconnected
918 * components.
919 *
920 * </para><para>
921 * Somewhat arbitrarily, igraph does not consider components containing
922 * a single vertex only as being biconnected. Isolated vertices will
923 * not be part of any of the biconnected components.
924 *
925 * \param graph The input graph.
926 * \param no The number of biconnected components will be stored here.
927 * \param tree_edges If not a NULL pointer, then the found components
928 * are stored here, in a list of vectors. Every vector in the list
929 * is a biconnected component, represented by its edges. More precisely,
930 * a spanning tree of the biconnected component is returned.
931 * Note you'll have to
932 * destroy each vector first by calling \ref igraph_vector_destroy()
933 * and then \ref igraph_free() on it, plus you need to call
934 * \ref igraph_vector_ptr_destroy() on the list to regain all
935 * allocated memory.
936 * \param component_edges If not a NULL pointer, then the edges of the
937 * biconnected components are stored here, in the same form as for
938 * \c tree_edges.
939 * \param components If not a NULL pointer, then the vertices of the
940 * biconnected components are stored here, in the same format as
941 * for the previous two arguments.
942 * \param articulation_points If not a NULL pointer, then the
943 * articulation points of the graph are stored in this vector.
944 * A vertex is an articulation point if its removal increases the
945 * number of (weakly) connected components in the graph.
946 * \return Error code.
947 *
948 * Time complexity: O(|V|+|E|), linear in the number of vertices and
949 * edges, but only if you do not calculate \c components and
950 * \c component_edges. If you calculate \c components, then it is
951 * quadratic in the number of vertices. If you calculate \c
952 * component_edges as well, then it is cubic in the number of
953 * vertices.
954 *
955 * \sa \ref igraph_articulation_points(), \ref igraph_clusters().
956 *
957 * \example examples/simple/igraph_biconnected_components.c
958 */
959
igraph_biconnected_components(const igraph_t * graph,igraph_integer_t * no,igraph_vector_ptr_t * tree_edges,igraph_vector_ptr_t * component_edges,igraph_vector_ptr_t * components,igraph_vector_t * articulation_points)960 int igraph_biconnected_components(const igraph_t *graph,
961 igraph_integer_t *no,
962 igraph_vector_ptr_t *tree_edges,
963 igraph_vector_ptr_t *component_edges,
964 igraph_vector_ptr_t *components,
965 igraph_vector_t *articulation_points) {
966
967 long int no_of_nodes = igraph_vcount(graph);
968 igraph_vector_long_t nextptr;
969 igraph_vector_long_t num, low;
970 igraph_vector_bool_t found;
971 igraph_vector_int_t *adjedges;
972 igraph_stack_t path;
973 igraph_vector_t edgestack;
974 igraph_inclist_t inclist;
975 long int i, counter, rootdfs = 0;
976 igraph_vector_long_t vertex_added;
977 long int comps = 0;
978 igraph_vector_ptr_t *mycomponents = components, vcomponents;
979
980 IGRAPH_CHECK(igraph_vector_long_init(&nextptr, no_of_nodes));
981 IGRAPH_FINALLY(igraph_vector_long_destroy, &nextptr);
982 IGRAPH_CHECK(igraph_vector_long_init(&num, no_of_nodes));
983 IGRAPH_FINALLY(igraph_vector_long_destroy, &num);
984 IGRAPH_CHECK(igraph_vector_long_init(&low, no_of_nodes));
985 IGRAPH_FINALLY(igraph_vector_long_destroy, &low);
986 IGRAPH_CHECK(igraph_vector_bool_init(&found, no_of_nodes));
987 IGRAPH_FINALLY(igraph_vector_bool_destroy, &found);
988
989 IGRAPH_CHECK(igraph_stack_init(&path, 100));
990 IGRAPH_FINALLY(igraph_stack_destroy, &path);
991 IGRAPH_VECTOR_INIT_FINALLY(&edgestack, 0);
992 IGRAPH_CHECK(igraph_vector_reserve(&edgestack, 100));
993
994 IGRAPH_CHECK(igraph_inclist_init(graph, &inclist, IGRAPH_ALL, IGRAPH_LOOPS_TWICE));
995 IGRAPH_FINALLY(igraph_inclist_destroy, &inclist);
996
997 IGRAPH_CHECK(igraph_vector_long_init(&vertex_added, no_of_nodes));
998 IGRAPH_FINALLY(igraph_vector_long_destroy, &vertex_added);
999
1000 if (no) {
1001 *no = 0;
1002 }
1003 if (tree_edges) {
1004 igraph_vector_ptr_clear(tree_edges);
1005 }
1006 if (components) {
1007 igraph_vector_ptr_clear(components);
1008 }
1009 if (component_edges) {
1010 igraph_vector_ptr_clear(component_edges);
1011 }
1012 if (articulation_points) {
1013 igraph_vector_clear(articulation_points);
1014 }
1015 if (component_edges && !components) {
1016 mycomponents = &vcomponents;
1017 IGRAPH_CHECK(igraph_vector_ptr_init(mycomponents, 0));
1018 IGRAPH_FINALLY(igraph_i_free_vectorlist, mycomponents);
1019 }
1020
1021 for (i = 0; i < no_of_nodes; i++) {
1022
1023 if (VECTOR(low)[i] != 0) {
1024 continue; /* already visited */
1025 }
1026
1027 IGRAPH_ALLOW_INTERRUPTION();
1028
1029 IGRAPH_CHECK(igraph_stack_push(&path, i));
1030 counter = 1;
1031 rootdfs = 0;
1032 VECTOR(low)[i] = VECTOR(num)[i] = counter++;
1033 while (!igraph_stack_empty(&path)) {
1034 long int n;
1035 long int act = (long int) igraph_stack_top(&path);
1036 long int actnext = VECTOR(nextptr)[act];
1037
1038 adjedges = igraph_inclist_get(&inclist, act);
1039 n = igraph_vector_int_size(adjedges);
1040 if (actnext < n) {
1041 /* Step down (maybe) */
1042 long int edge = (long int) VECTOR(*adjedges)[actnext];
1043 long int nei = IGRAPH_OTHER(graph, edge, act);
1044 if (VECTOR(low)[nei] == 0) {
1045 if (act == i) {
1046 rootdfs++;
1047 }
1048 IGRAPH_CHECK(igraph_vector_push_back(&edgestack, edge));
1049 IGRAPH_CHECK(igraph_stack_push(&path, nei));
1050 VECTOR(low)[nei] = VECTOR(num)[nei] = counter++;
1051 } else {
1052 /* Update low value if needed */
1053 if (VECTOR(num)[nei] < VECTOR(low)[act]) {
1054 VECTOR(low)[act] = VECTOR(num)[nei];
1055 }
1056 }
1057 VECTOR(nextptr)[act] += 1;
1058 } else {
1059 /* Step up */
1060 igraph_stack_pop(&path);
1061 if (!igraph_stack_empty(&path)) {
1062 long int prev = (long int) igraph_stack_top(&path);
1063 /* Update LOW value if needed */
1064 if (VECTOR(low)[act] < VECTOR(low)[prev]) {
1065 VECTOR(low)[prev] = VECTOR(low)[act];
1066 }
1067 /* Check for articulation point */
1068 if (VECTOR(low)[act] >= VECTOR(num)[prev]) {
1069 if (articulation_points && !VECTOR(found)[prev]
1070 && prev != i /* the root */) {
1071 IGRAPH_CHECK(igraph_vector_push_back(articulation_points, prev));
1072 VECTOR(found)[prev] = 1;
1073 }
1074 if (no) {
1075 *no += 1;
1076 }
1077
1078 /*------------------------------------*/
1079 /* Record the biconnected component just found */
1080 if (tree_edges || mycomponents) {
1081 igraph_vector_t *v = 0, *v2 = 0;
1082 comps++;
1083 if (tree_edges) {
1084 v = IGRAPH_CALLOC(1, igraph_vector_t);
1085 if (!v) {
1086 IGRAPH_ERROR("Out of memory", IGRAPH_ENOMEM);
1087 }
1088 IGRAPH_CHECK(igraph_vector_init(v, 0));
1089 IGRAPH_FINALLY(igraph_vector_destroy, v);
1090 }
1091 if (mycomponents) {
1092 v2 = IGRAPH_CALLOC(1, igraph_vector_t);
1093 if (!v2) {
1094 IGRAPH_ERROR("Out of memory", IGRAPH_ENOMEM);
1095 }
1096 IGRAPH_CHECK(igraph_vector_init(v2, 0));
1097 IGRAPH_FINALLY(igraph_vector_destroy, v2);
1098 }
1099
1100 while (!igraph_vector_empty(&edgestack)) {
1101 long int e = (long int) igraph_vector_pop_back(&edgestack);
1102 long int from = IGRAPH_FROM(graph, e);
1103 long int to = IGRAPH_TO(graph, e);
1104 if (tree_edges) {
1105 IGRAPH_CHECK(igraph_vector_push_back(v, e));
1106 }
1107 if (mycomponents) {
1108 if (VECTOR(vertex_added)[from] != comps) {
1109 VECTOR(vertex_added)[from] = comps;
1110 IGRAPH_CHECK(igraph_vector_push_back(v2, from));
1111 }
1112 if (VECTOR(vertex_added)[to] != comps) {
1113 VECTOR(vertex_added)[to] = comps;
1114 IGRAPH_CHECK(igraph_vector_push_back(v2, to));
1115 }
1116 }
1117 if (from == prev || to == prev) {
1118 break;
1119 }
1120 }
1121
1122 if (mycomponents) {
1123 IGRAPH_CHECK(igraph_vector_ptr_push_back(mycomponents, v2));
1124 IGRAPH_FINALLY_CLEAN(1);
1125 }
1126 if (tree_edges) {
1127 IGRAPH_CHECK(igraph_vector_ptr_push_back(tree_edges, v));
1128 IGRAPH_FINALLY_CLEAN(1);
1129 }
1130 if (component_edges) {
1131 igraph_vector_t *nodes = VECTOR(*mycomponents)[comps - 1];
1132 igraph_vector_t *vv = IGRAPH_CALLOC(1, igraph_vector_t);
1133 long int ii, no_vert = igraph_vector_size(nodes);
1134 if (!vv) {
1135 IGRAPH_ERROR("Out of memory", IGRAPH_ENOMEM);
1136 }
1137 IGRAPH_CHECK(igraph_vector_init(vv, 0));
1138 IGRAPH_FINALLY(igraph_vector_destroy, vv);
1139 for (ii = 0; ii < no_vert; ii++) {
1140 long int vert = (long int) VECTOR(*nodes)[ii];
1141 igraph_vector_int_t *edges = igraph_inclist_get(&inclist,
1142 vert);
1143 long int j, nn = igraph_vector_int_size(edges);
1144 for (j = 0; j < nn; j++) {
1145 long int e = (long int) VECTOR(*edges)[j];
1146 long int nei = IGRAPH_OTHER(graph, e, vert);
1147 if (VECTOR(vertex_added)[nei] == comps && nei < vert) {
1148 IGRAPH_CHECK(igraph_vector_push_back(vv, e));
1149 }
1150 }
1151 }
1152 IGRAPH_CHECK(igraph_vector_ptr_push_back(component_edges, vv));
1153 IGRAPH_FINALLY_CLEAN(1);
1154 }
1155 } /* record component if requested */
1156 /*------------------------------------*/
1157
1158 }
1159 } /* !igraph_stack_empty(&path) */
1160 }
1161
1162 } /* !igraph_stack_empty(&path) */
1163
1164 if (articulation_points && rootdfs >= 2) {
1165 IGRAPH_CHECK(igraph_vector_push_back(articulation_points, i));
1166 }
1167
1168 } /* i < no_of_nodes */
1169
1170 if (mycomponents != components) {
1171 igraph_i_free_vectorlist(mycomponents);
1172 IGRAPH_FINALLY_CLEAN(1);
1173 }
1174
1175 igraph_vector_long_destroy(&vertex_added);
1176 igraph_inclist_destroy(&inclist);
1177 igraph_vector_destroy(&edgestack);
1178 igraph_stack_destroy(&path);
1179 igraph_vector_bool_destroy(&found);
1180 igraph_vector_long_destroy(&low);
1181 igraph_vector_long_destroy(&num);
1182 igraph_vector_long_destroy(&nextptr);
1183 IGRAPH_FINALLY_CLEAN(8);
1184
1185 return 0;
1186 }
1187
1188
1189 /* igraph_bridges -- find all bridges in the graph */
1190 /* The algorithm is based on https://www.geeksforgeeks.org/bridge-in-a-graph/
1191 but instead of keeping track of the parent of each vertex in the DFS tree
1192 we keep track of its incoming edge. This is necessary to support multigraphs. */
1193
igraph_i_bridges_rec(const igraph_t * graph,const igraph_inclist_t * il,igraph_integer_t u,igraph_integer_t * time,igraph_vector_t * bridges,igraph_vector_bool_t * visited,igraph_vector_int_t * disc,igraph_vector_int_t * low,igraph_vector_int_t * incoming_edge)1194 static int igraph_i_bridges_rec(
1195 const igraph_t *graph, const igraph_inclist_t *il, igraph_integer_t u,
1196 igraph_integer_t *time, igraph_vector_t *bridges, igraph_vector_bool_t *visited,
1197 igraph_vector_int_t *disc, igraph_vector_int_t *low, igraph_vector_int_t *incoming_edge)
1198 {
1199 igraph_vector_int_t *incedges;
1200 long nc; /* neighbour count */
1201 long i;
1202
1203 VECTOR(*visited)[u] = 1;
1204
1205 *time += 1;
1206
1207 VECTOR(*disc)[u] = *time;
1208 VECTOR(*low)[u] = *time;
1209
1210 incedges = igraph_inclist_get(il, u);
1211 nc = igraph_vector_int_size(incedges);
1212 for (i = 0; i < nc; ++i) {
1213 long edge = (long) VECTOR(*incedges)[i];
1214 igraph_integer_t v = IGRAPH_TO(graph, edge) == u ? IGRAPH_FROM(graph, edge) : IGRAPH_TO(graph, edge);
1215
1216 if (! VECTOR(*visited)[v]) {
1217 VECTOR(*incoming_edge)[v] = edge;
1218 IGRAPH_CHECK(igraph_i_bridges_rec(graph, il, v, time, bridges, visited, disc, low, incoming_edge));
1219
1220 VECTOR(*low)[u] = VECTOR(*low)[u] < VECTOR(*low)[v] ? VECTOR(*low)[u] : VECTOR(*low)[v];
1221
1222 if (VECTOR(*low)[v] > VECTOR(*disc)[u]) {
1223 IGRAPH_CHECK(igraph_vector_push_back(bridges, edge));
1224 }
1225 } else if (edge != VECTOR(*incoming_edge)[u]) {
1226 VECTOR(*low)[u] = VECTOR(*low)[u] < VECTOR(*disc)[v] ? VECTOR(*low)[u] : VECTOR(*disc)[v];
1227 }
1228 }
1229
1230 return IGRAPH_SUCCESS;
1231 }
1232
1233 /**
1234 * \function igraph_bridges
1235 * Find all bridges in a graph.
1236 *
1237 * An edge is a bridge if its removal increases the number of (weakly)
1238 * connected components in the graph.
1239 *
1240 * \param graph The input graph.
1241 * \param res Pointer to an initialized vector, the
1242 * bridges will be stored here as edge indices.
1243 * \return Error code.
1244 *
1245 * Time complexity: O(|V|+|E|), linear in the number of vertices and edges.
1246 *
1247 * \sa \ref igraph_articulation_points(), \ref igraph_biconnected_components(), \ref igraph_clusters()
1248 */
1249
igraph_bridges(const igraph_t * graph,igraph_vector_t * bridges)1250 int igraph_bridges(const igraph_t *graph, igraph_vector_t *bridges) {
1251 igraph_inclist_t il;
1252 igraph_vector_bool_t visited;
1253 igraph_vector_int_t disc, low;
1254 igraph_vector_int_t incoming_edge;
1255 long n;
1256 long i;
1257 igraph_integer_t time;
1258
1259 n = igraph_vcount(graph);
1260
1261 IGRAPH_CHECK(igraph_inclist_init(graph, &il, IGRAPH_ALL, IGRAPH_LOOPS_TWICE));
1262 IGRAPH_FINALLY(igraph_inclist_destroy, &il);
1263
1264 IGRAPH_CHECK(igraph_vector_bool_init(&visited, n));
1265 IGRAPH_FINALLY(igraph_vector_bool_destroy, &visited);
1266
1267 IGRAPH_CHECK(igraph_vector_int_init(&disc, n));
1268 IGRAPH_FINALLY(igraph_vector_int_destroy, &disc);
1269
1270 IGRAPH_CHECK(igraph_vector_int_init(&low, n));
1271 IGRAPH_FINALLY(igraph_vector_int_destroy, &low);
1272
1273 IGRAPH_CHECK(igraph_vector_int_init(&incoming_edge, n));
1274 IGRAPH_FINALLY(igraph_vector_int_destroy, &incoming_edge);
1275 for (i = 0; i < n; ++i) {
1276 VECTOR(incoming_edge)[i] = -1;
1277 }
1278
1279 igraph_vector_clear(bridges);
1280
1281 time = 0;
1282 for (i = 0; i < n; ++i)
1283 if (! VECTOR(visited)[i]) {
1284 IGRAPH_CHECK(igraph_i_bridges_rec(graph, &il, i, &time, bridges, &visited, &disc, &low, &incoming_edge));
1285 }
1286
1287 igraph_vector_int_destroy(&incoming_edge);
1288 igraph_vector_int_destroy(&low);
1289 igraph_vector_int_destroy(&disc);
1290 igraph_vector_bool_destroy(&visited);
1291 igraph_inclist_destroy(&il);
1292 IGRAPH_FINALLY_CLEAN(5);
1293
1294 return IGRAPH_SUCCESS;
1295 }
1296
1297 /**
1298 * \ingroup structural
1299 * \function igraph_subcomponent
1300 * \brief The vertices in the same component as a given vertex.
1301 *
1302 * \param graph The graph object.
1303 * \param res The result, vector with the ids of the vertices in the
1304 * same component.
1305 * \param vertex The id of the vertex of which the component is
1306 * searched.
1307 * \param mode Type of the component for directed graphs, possible
1308 * values:
1309 * \clist
1310 * \cli IGRAPH_OUT
1311 * the set of vertices reachable \em from the
1312 * \p vertex,
1313 * \cli IGRAPH_IN
1314 * the set of vertices from which the
1315 * \p vertex is reachable.
1316 * \cli IGRAPH_ALL
1317 * the graph is considered as an
1318 * undirected graph. Note that this is \em not the same
1319 * as the union of the previous two.
1320 * \endclist
1321 * \return Error code:
1322 * \clist
1323 * \cli IGRAPH_ENOMEM
1324 * not enough memory for temporary data.
1325 * \cli IGRAPH_EINVVID
1326 * \p vertex is an invalid vertex id
1327 * \cli IGRAPH_EINVMODE
1328 * invalid mode argument passed.
1329 * \endclist
1330 *
1331 * Time complexity: O(|V|+|E|),
1332 * |V| and
1333 * |E| are the number of vertices and
1334 * edges in the graph.
1335 *
1336 * \sa \ref igraph_induced_subgraph() if you want a graph object consisting only
1337 * a given set of vertices and the edges between them.
1338 */
igraph_subcomponent(const igraph_t * graph,igraph_vector_t * res,igraph_real_t vertex,igraph_neimode_t mode)1339 int igraph_subcomponent(const igraph_t *graph, igraph_vector_t *res, igraph_real_t vertex,
1340 igraph_neimode_t mode) {
1341
1342 long int no_of_nodes = igraph_vcount(graph);
1343 igraph_dqueue_t q = IGRAPH_DQUEUE_NULL;
1344 char *already_added;
1345 long int i, vsize;
1346 igraph_vector_t tmp = IGRAPH_VECTOR_NULL;
1347
1348 if (!IGRAPH_FINITE(vertex) || vertex < 0 || vertex >= no_of_nodes) {
1349 IGRAPH_ERROR("Vertex id out of range.", IGRAPH_EINVVID);
1350 }
1351 if (mode != IGRAPH_OUT && mode != IGRAPH_IN &&
1352 mode != IGRAPH_ALL) {
1353 IGRAPH_ERROR("Invalid mode argument.", IGRAPH_EINVMODE);
1354 }
1355
1356 already_added = IGRAPH_CALLOC(no_of_nodes, char);
1357 if (already_added == 0) {
1358 IGRAPH_ERROR("Subcomponent failed.", IGRAPH_ENOMEM);
1359 }
1360 IGRAPH_FINALLY(igraph_free, already_added);
1361
1362 igraph_vector_clear(res);
1363
1364 IGRAPH_VECTOR_INIT_FINALLY(&tmp, 0);
1365 IGRAPH_DQUEUE_INIT_FINALLY(&q, 100);
1366
1367 IGRAPH_CHECK(igraph_dqueue_push(&q, vertex));
1368 IGRAPH_CHECK(igraph_vector_push_back(res, vertex));
1369 already_added[(long int)vertex] = 1;
1370
1371 while (!igraph_dqueue_empty(&q)) {
1372 long int actnode = (long int) igraph_dqueue_pop(&q);
1373
1374 IGRAPH_ALLOW_INTERRUPTION();
1375
1376 IGRAPH_CHECK(igraph_neighbors(graph, &tmp, (igraph_integer_t) actnode,
1377 mode));
1378 vsize = igraph_vector_size(&tmp);
1379 for (i = 0; i < vsize; i++) {
1380 long int neighbor = (long int) VECTOR(tmp)[i];
1381
1382 if (already_added[neighbor]) {
1383 continue;
1384 }
1385 already_added[neighbor] = 1;
1386 IGRAPH_CHECK(igraph_vector_push_back(res, neighbor));
1387 IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
1388 }
1389 }
1390
1391 igraph_dqueue_destroy(&q);
1392 igraph_vector_destroy(&tmp);
1393 IGRAPH_FREE(already_added);
1394 IGRAPH_FINALLY_CLEAN(3);
1395
1396 return IGRAPH_SUCCESS;
1397 }
1398