summaryrefslogtreecommitdiff
path: root/src/johnson/johnson_cycles.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/johnson/johnson_cycles.c')
-rw-r--r--src/johnson/johnson_cycles.c328
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);
+}
generated by cgit v1.2.3 (git 2.25.1) at 2024-05-04 05:30:54 +0000