diff options
Diffstat (limited to 'src/johnson/johnson_cycles.c')
-rw-r--r-- | src/johnson/johnson_cycles.c | 328 |
1 files changed, 328 insertions, 0 deletions
diff --git a/src/johnson/johnson_cycles.c b/src/johnson/johnson_cycles.c new file mode 100644 index 0000000..5135b20 --- /dev/null +++ b/src/johnson/johnson_cycles.c @@ -0,0 +1,328 @@ +/** + * This program is free software: you can redistribute it and/or + * modify it under the terms of the GNU General Public License as + * published by the Free Software Foundation, either version 3 of the + * License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see + * <http://www.gnu.org/licenses/>. + * + * (c) Vincenzo Nicosia 2009-2017 -- <v.nicosia@qmul.ac.uk> + * + * This file is part of NetBunch, a package for complex network + * analysis and modelling. For more information please visit: + * + * http://www.complex-networks.net/ + * + * If you use this software, please add a reference to + * + * V. Latora, V. Nicosia, G. Russo + * "Complex Networks: Principles, Methods and Applications" + * Cambridge University Press (2017) + * ISBN: 9781107103184 + * + *********************************************************************** + * + * Compute all the simple cycles of an undirected graph provided as + * input (optionally, only up to a certain length) using the Johnson + * algorithm. + * + * References: + * + * [1] D. B. Johnson. "Finding All the Elementary Circuits of a + * Directed Graph". SIAM J. Comput. 4 (1975), 77-84. + * + */ + +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "utils.h" +#include "stack.h" + + +#define FALSE 0 +#define TRUE 1 + + +void usage(char *argv[]){ + + printf("********************************************************************\n" + "** **\n" + "** -*- johnson_cycles -*- **\n" + "** **\n" + "** Enumerates all the elementary cycles of a graph on input. **\n" + "** **\n" + "** If 'graph_in' is equal to '-' (dash), read the edge list **\n" + "** from standard input (STDIN). **\n" + "** **\n" + "** If 'max_length' is provided, ignore all cycles whose length **\n" + "** is larger than 'max_length' **\n" + "** **\n" + "** The program prints on output a row for each cycle length, **\n" + "** in the format: **\n" + "** **\n" + "** 2 N_2 **\n" + "** 3 N_3 **\n" + "** 4 N_4 **\n" + "** 5 N_5 **\n" + "** .... **\n" + "** **\n" + "** where 2, 3, 4, 5, ... is the length of the cycle, and N_2, **\n" + "** N_3, N_4, N_5, ... is the number of cycles of that length. **\n" + "** Notice that in undirected graphs each cycle is reported **\n" + "** (and counted) twice, once for each possible direction in **\n" + "** which the cycle can be traversed. **\n" + "** **\n" + "** If the third parameter is equal to 'SHOW', the program will **\n" + "** dump all the found cycles on STDERR, using the format: **\n" + "** **\n" + "** N_(l-1) N_(l-2)... N_0 **\n" + "** **\n" + "** where N_(l-i) is the i-th node in the cycle starting at node **\n" + "** N_0. Notice that according to this notation we report only **\n" + "** 'l' nodes belonging to the cycle, in reversed traversal **\n" + "** order, avoiding to report the starting node (N_0) twice. **\n" + "** **\n" + "** The program shows also cycles of length '2', that in **\n" + "** undirected graphs correspond to the degenerate cycle obtained **\n" + "** by traversing the same edge in the two directions. **\n" + "** **\n" + "********************************************************************\n" + " This is Free Software - You can use and distribute it under \n" + " the terms of the GNU General Public License, version 3 or later\n\n" + " Please visit http://www.complex-networks.net for more information\n\n" + " (c) Vincenzo Nicosia 2009-2017 (v.nicosia@qmul.ac.uk)\n" + "********************************************************************\n\n" + ); + printf("Usage: %s <graph_in> [<max_length> [SHOW]]\n", argv[0]); + +} + + +typedef struct{ + unsigned int num; + double *c; +} count_t; + + +/** + * + * This structure contains all the information needed by the algorithm + * to run: + * + * - N is the number of nodes + * + * - K is twice the number of edges + * + * - J_slap and r_slap are the usual SLAP representation of the + * adjacency matrix + * + * - B_len is an array with the length of the set B(u) for each node u + * + * - B is a vector of character of the same size of J_slap, which + * stores the values of the indicator vector B(u) for each node u + * + * - blocked is an indicator vector of size N + * + * - v is the stack of nodes + * + * - stats maintains the counts for all the cycle lengths + * + */ +typedef struct{ + unsigned int N; + unsigned int K; + unsigned int *J_slap; + unsigned int *r_slap; + my_stack_t *stack; + int max_depth; + count_t stats; + unsigned int *B_len; + unsigned int *B; + char *blocked; + char show; +} algo_info_t; + + + + +void unblock(algo_info_t *info, unsigned int u){ + + unsigned int w, idx; + + info->blocked[u] = FALSE; + + while(info->B_len[u] >0){ + info->B_len[u] -= 1; + idx = info->r_slap[u] + info->B_len[u]; + w = info->B[idx]; + if (info->blocked[w] == TRUE){ + unblock(info, w); + } + } +} + + + + +/* add v to B(w) */ +void set_B(algo_info_t *info, unsigned int v, unsigned int w){ + + unsigned int idx, i; + + idx = info->r_slap[w]; + + + for(i=0; i <info->B_len[w]; i ++){ + if (info->B[idx + i] == v) + return; + } + info->B[idx + info->B_len[w]] = v; + info->B_len[w] += 1; + +} + +char circuit(algo_info_t *info, unsigned int s, unsigned int v){ + + unsigned int i, w; + char f = FALSE; + + if (v < s) + return FALSE; + + /* check if maximum depth has been reached */ + if(stack_size(info->stack) == info->max_depth){ + return TRUE; + } + + stack_push(info->stack, v); + + info->blocked[v] = TRUE; + for(i=info->r_slap[v]; i< info->r_slap[v+1]; i++){ + w = info->J_slap[i]; + if (w == s){ + /* output a cycle here */ + if (info->show) + stack_dump(info->stack, stderr); + info->stats.c[stack_size(info->stack)] += 1; + f = TRUE; + } + if (w < s) + continue; /* We only consider nodes with higher IDs */ + else + if (! (info->blocked[w])){ + if (circuit(info, s, w) == TRUE) + f = TRUE; + } + } + if (f == TRUE) { + unblock(info, v); + } + else{ + for(i=info->r_slap[v]; i < info->r_slap[v+1]; i++){ + w = info->J_slap[i]; + if (w > s) + set_B(info, v, w); /* add v to B(w) only if v is not already in B(w) */ + } + } + stack_pop(info->stack, &w); + return f; +} + + +void johnson_cycles(algo_info_t *info){ + + unsigned int i, N; + + N = info->N; + + + for(i=0; i<N; i++){ + /* clear the vector of blocked nodes */ + memset(info->blocked, 0, N * sizeof(char)); + + /* clear the indicator vector B */ + + /* we reset the length of each set B*/ + memset(info->B_len, 0, N * sizeof(unsigned int)); + circuit(info, i, i); + } + + +} + + +int main(int argc, char *argv[]){ + + algo_info_t info; + FILE *filein; + unsigned int i; + + if (argc < 2){ + usage(argv); + exit(1); + } + + info.J_slap = NULL; + info.r_slap = NULL; + + + if(!strcmp(argv[1], "-")){ + filein = stdin; + } + else{ + filein = openfile_or_exit(argv[1], "r", 2); + } + + read_slap(filein, &(info.K), &(info.N), &(info.J_slap), &(info.r_slap)); + + fclose(filein); + + info.B = malloc(info.r_slap[info.N] * sizeof(unsigned int)); /* This is the indicator vector B(u) */ + info.B_len = malloc(info.N * sizeof(unsigned int)); /* lengths of the B(u) */ + info.blocked = malloc(info.N * sizeof(char)); /* This is the indicator vector + of blocked nodes */ + + info.stack = malloc(sizeof(my_stack_t)); + stack_create(info.stack); + + info.show = FALSE; + if (argc > 2){ + info.max_depth = atoi(argv[2]); + if (argc > 3){ + if (!my_strcasecmp(argv[3], "show")){ + info.show = TRUE; + } + } + } + else{ + info.max_depth = info.N; + } + + /* Here we initialise the counters */ + info.stats.num = info.max_depth + 1; + info.stats.c = malloc(info.stats.num * sizeof(double)); + /* Here we set the counters to zero */ + memset(info.stats.c, 0, info.stats.num * sizeof(double)); + johnson_cycles(&info); + + for(i=2; i<info.stats.num ; i ++){ + printf("%d %1.20g\n", i, info.stats.c[i]); + } + free(info.B); + free(info.B_len); + free(info.blocked); + free(info.stack->v); + free(info.stack); + free(info.stats.c); + free(info.r_slap); + free(info.J_slap); +} |