Added function to define of aglomerate type

statgen/Makefile

@@ -3,7 +3,7 @@ PROJECT=STATGEN
 CC=gcc
 CFLAGS=-c -Wall -fPIC
 LDFLAGS=-lm
-SOURCES=main.c coords.c int2str.c messages.c stat_print.c stat_select.c stat_sort.c summary_stat.c
+SOURCES=main.c coords.c graph.c int2str.c messages.c stat_print.c stat_select.c stat_sort.c summary_stat.c
 OBJECTS=$(SOURCES:.c=.o)
 EXECUTABLE=statgen
 

statgen/graph.c

@@ -0,0 +1,144 @@
+/* Library for graph structure analyze
+ * Usage:
+ *                                         graph_analyze (N, connect, max_depth)
+ */
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+
+int check_cycle (const int N, const int *matrix)
+// function to return number of cycles
+{
+  int cycle, i, j;
+/* cycle            - number of cycle
+ */
+  
+  cycle = 0;
+  for (i=0; i<N; i++)
+    for (j=0; j<N; j++)
+      cycle += matrix[i*N+j];
+  
+//   for linear (0.5*cycle == N-1)
+  cycle = 0.5 * cycle - (N - 1);
+  
+  return cycle;
+}
+
+
+int check_cycle_size (const int N, const int *matrix, const int size)
+// function to return number of cycles of certain size
+{
+  int cur_N, cycle, i, j, k, n, num, *submatrix, *vertex;
+/* cur_N            - current number of elements in submatrix
+ * cycle            - if ((0.5 * cycle) == size) that cycle exist
+ * n                - number of samples
+ * num              - number of cycles of certain size
+ * submatrix        - connectivity matrix for subgraph
+ * vertex           - vertexes of subgraph
+ */
+  
+  num = 0;
+  submatrix = (int *) malloc (size * size * sizeof (int));
+  vertex = (int *) malloc (N * sizeof (int));
+  
+//   matrix generation from
+//   http://wincode.org/acm-icpc/subsets-generation
+  n = pow (2, N);
+  for (i=0; i<n; i++)
+  {
+    cur_N = 0;
+    for (j=0; j<N; j++)
+      if ( i & (1 << j))
+      {
+        vertex[cur_N] = j;
+        cur_N++;
+      }
+    
+    if (cur_N == size)
+    {
+//       copy connectivity matrix
+      for (j=0; j<cur_N; j++)
+        for (k=0; k<cur_N; k++)
+          submatrix[j*cur_N+k] = matrix[vertex[j]*N+vertex[k]];
+      
+      cycle = 0;
+//       analyze subgraph
+      for (j=0; j<cur_N; j++)
+        for (k=0; k<cur_N; k++)
+          cycle += submatrix[j*cur_N+j];
+      
+      if ((0.5 * cycle) == size)
+        num++;
+    }
+  }
+  
+  free (vertex);
+  free (submatrix);
+  
+  return num;
+}
+
+
+int check_tail (const int N, const int *matrix)
+// function to return number of tails
+{
+  int i, j, pn, tail;
+/* pn               - weight one vertex
+ * tail             - number of tails
+ */
+  
+  tail = 0;
+  for (i=0; i<N; i++)
+  {
+    pn = 0;
+    for (j=0; j<N; j++)
+      pn += matrix[i*N+j];
+    if (pn == 1)
+      tail++;
+  }
+  
+  return tail;
+}
+
+
+int graph_analyze (const int N, const int *matrix, const int max_depth, char *iso)
+/* N                - number of vertex in graph
+ * matrix           - connectivity matrix
+ * max_depth        - maximum depth for check_cycle_size
+ * iso              - isomorphism class
+ */
+{
+  int cycle, depth, i, *n_cycle, tail;
+/* cycle            - number of cycles
+ * depth            - depth for check_cycle_size
+ * n_cycle          - number of cycles of certain size
+ * tail             - number of tails
+ */
+  
+  if (max_depth > N)
+    depth = N;
+  else
+    depth = max_depth;
+  n_cycle = (int *) malloc ((max_depth-2) * sizeof (int));
+  
+  tail = check_tail (N, matrix);
+  cycle = check_cycle (N, matrix);
+  if (cycle > 0)
+    for (i=0; i<depth-2; i++)
+      n_cycle[i] = check_cycle_size (N, matrix, i+3);
+  else
+    for (i=0; i<depth-2; i++)
+      n_cycle[i] = 0;
+  for (i=depth-2; i<max_depth-2; i++)
+    n_cycle[i] = 0;
+  
+  sprintf (iso, "%i.%i.%i", N, tail, cycle);
+  for (i=0; i<max_depth-2; i++)
+    sprintf (iso, "%s.%i", iso, n_cycle[i]);
+  
+  free (n_cycle);
+  
+  return 1;
+}

statgen/stat_print.c

@@ -10,6 +10,10 @@
 #include <stdlib.h>
 
 
+// prototype
+int graph_analyze (const int, const int *, const int, char *);
+
+
 int printing_agl (const char *input, const char *output, const int *connect, 
                   const int num_mol, const int *true_label_mol, const int *num_mol_agl, 
                   const int *agl, const int *stat, int *type_agl)
@@ -24,10 +28,10 @@ int printing_agl (const char *input, const char *output, const int *connect,
  * type_agl         - massive of numbers of aglomerate types
  */
 {
-  int i, iso, j, k, p, type, *label_matrix, **matrix;
+  char iso[256];
+  int i, j, k, *label_matrix, *matrix;
   FILE *f_out;
 /* iso              - isomorphic graph in database
- * type             - number of cycle in aglomerates
  * label_matrix     - massive of indexes of molecule
  * matrix           - connectivity graph
  * f_out            - output file
@@ -47,13 +51,10 @@ int printing_agl (const char *input, const char *output, const int *connect,
     if (num_mol_agl[i] > 0)
     {
 //       creating connectivity graph
-      matrix = (int **) malloc (num_mol_agl[i] * sizeof (int *));
+      matrix = (int *) malloc (num_mol_agl[i] * num_mol_agl[i] * sizeof (int));
       for (j=0; j<num_mol_agl[i]; j++)
-      {
-        matrix[j] = (int *) malloc (num_mol_agl[i] * sizeof (int));
         for (k=0; k<num_mol_agl[i]; k++)
-          matrix[j][k] = 0;
-      }
+          matrix[num_mol_agl[i]*j+k] = 0;
       label_matrix = (int *) malloc (num_mol * sizeof (int));
       if ((matrix == NULL) ||
         (label_matrix == NULL))
@@ -65,42 +66,27 @@ int printing_agl (const char *input, const char *output, const int *connect,
         for (k=j+1; k<num_mol_agl[i]; k++)
           if (connect[num_mol*agl[num_mol*i+j]+agl[num_mol*i+k]] == 1)
           {
-            matrix[label_matrix[agl[num_mol*i+j]]][label_matrix[agl[num_mol*i+k]]] = 1;
-            matrix[label_matrix[agl[num_mol*i+k]]][label_matrix[agl[num_mol*i+j]]] = 1;
-          }
-//       TODO: analyze of topology
-      iso = 0;
-      p = 0;
-      for (j=0; j<num_mol_agl[i]; j++)
-        for (k=0; k<num_mol_agl[i]; k++)
-          p += matrix[j][k];
-      if (p == (2*num_mol_agl[i]-2))
-      {
-        type = 0;
-        type_agl[0]++;
-      }
-      else
-      {
-        type = (p - (2*num_mol_agl[i]-2)) / 2;
-        type_agl[1]++;
+            matrix[label_matrix[agl[num_mol*i+j]]*num_mol_agl[i]+label_matrix[agl[num_mol*i+k]]] = 1;
+            matrix[label_matrix[agl[num_mol*i+k]]*num_mol_agl[i]+label_matrix[agl[num_mol*i+j]]] = 1;
           }
 
+//       graph topology analyze
+      graph_analyze (num_mol_agl[i], matrix, 3, iso);
+      
 //       printing class of aglomerate
-      fprintf (f_out, "AGL=%i=%i=%i\n", num_mol_agl[i], type, iso);
+      fprintf (f_out, "AGL=%i=%s\n", num_mol_agl[i], iso);
       for (j=0; j<num_mol_agl[i]; j++)
       {
         fprintf (f_out, "%7i=", true_label_mol[agl[num_mol*i+j]]);
         for (k=0; k<num_mol_agl[i]; k++)
         {
-          if (matrix[j][k] == 1)
+          if (matrix[j*num_mol_agl[i]+k] == 1)
             fprintf (f_out, "%i,", true_label_mol[agl[num_mol*i+k]]);
         }
         fprintf (f_out, "\n");
       }
       
 //       free memory
-      for (j=0; j<num_mol_agl[i]; j++)
-        free (matrix[j]);
       free (matrix);
       free (label_matrix);
     }