Linear average time automorphism algorithm for random graphs.

Master’s Paper · 1996

Luis Ferreira · Department of Physics, Pennsylvania State University

This is my 1996 Master’s paper, reproduced here in full and unedited. The algorithm computes the Variety of a random graph — a measure, in Lee Smolin’s sense, of how well each point of a structure can be told apart from every other. The three figures are new: they were recomputed in 2026 with Python and NetworkX to illustrate the text, and did not appear in the original document. The original C source is at the end of the post and on GitHub.

Abstract

An algorithm is developed and implemented, whose purpose is to determine the automorphisms of a random graph in linear average time and to calculate a related function of a graph called the Variety. For most part of the probability range the Algorithm is very efficient and so we are able to compute, in a couple of minutes, graphs with up to 1000 vertices. For a short range of probabilities close to zero (or on the symmetric position in relation to 0.5) it might take a much longer time, due, probably, to an increase in symmetry and disconnectivity of the graph. This range's width depends inversely on the number of vertices of the graph.

Problem

There is no know solution for the isomorphism problem for all graphs. It might be a NP-complete problem, meaning not only that it might not be solvable in polynomial time but also that it might belong to a class of problems which if one is proven solvable then all of them are [1, 2]. However there are strategies that can be implemented to solve special cases like the random graph automorphism.

Physical Motivation

The concept of variety, as defined in the following discussion, may be useful for the study of complexity in diffeomorphism invariant theories like gravitational and cosmological theories. It may also provide a quantitative tool to study self-organizing systems [3, 4, 5].

Introduction

A graph is a set of vertices and connections between pairs of vertices. A connection can be in either one of two states: on or off. In a random graph each connection has a probability p of being on (p is between 0 and 1).

The 1st order view of a vertex is the set of vertices to which it is connected.

The nth order view of a vertex is the sub-graph formed by all vertices within n steps from the original vertex and all the connections among them.

The Indifference of a vertex is the number of steps we must travel to obtain a view that is different from all the other views. If we cannot distinguish a vertex from another we say that they are isomorphic and then they have Indifference equal to infinite.

The Indifference of a graph is the sum of the indifference of all the vertices.

We want to calculate a function of a random graph called the Variety.

The Variety measures the difference between the view of each point and all the other points in the graph. It is the normalised sum of the inverses of all the Indifference.

Heuristically this Variety is a measure of "the ability to differentiate each point of the graph from all the others".

(The definitions of this section are in accordance with [3, 4, 5]).

Random graph with vertices coloured by indifference depth
Figure 1 — Indifference across a random graph. G(N = 120, p = 0.04). Each vertex is coloured by its indifference: the number of refinement steps needed before its view becomes unique. Bright vertices are distinguished almost immediately; dark ones sit in symmetric pockets of the graph and resist separation for several steps.

Preliminaries

The program requires some input, namely the probability p that each connection is on and the number of runs (ensemble of random graphs).

The resulting Variety data is written in a file called "var.out" and on the screen as well.

It creates a random graph with probability p and stores it as an Adjacency Matrix.

An Adjacency Matrix is a matrix of zeros and ones whose rows and columns represent the vertices and whose entries represent the connections.

The program accepts a probability range between 0 and 1 but it treats the graphs symmetrically on each side of the value 0.5 (the variety will be the same if we exchange the 0 and 1 entries of the adjacency matrix).

The way the program works is by dividing iteratively the vertices among different baskets which correspond to different views.

Each iteration extracts the n th order view of each vertex (sub-graph X).

Once the program obtains a basket with a single vertex it calculates the variety of that vertex and adds it to the total variety. It proceeds until it has calculated the variety of the entire graph.

To test for isomorphism (same view) of different vertices we adapted an algorithm by Babai and Kucera [6].

The algorithm attributes a Canonical Label to each subgraph.

A Canonical Label is an assignment of a linear order to the set of vertices of a graph.

Two graphs are isomorphic if there exists an order preserving map between their respective Canonical Labels.

Canonical Label Algorithm

(in accordance with [6]).

A single vertex classification procedure (procedure A) has a very high probability of success after just two refinement steps (since we are dealing with random graphs).

An additional depth-first search (procedure B) yields a canonical labelling of all graphs.

V is the set of vertices.

A partial order is a an order of equivalence classes of vertices, C1 < C2 < ... < Cn.

A Canonical Vertex Classification is a partial order P.

Procedure A: Vertex classification procedure.

The input is a sub-graph X.

The first iteration of the classification procedure orders the vertices into equivalence classes according to their valences.

The second iteration refines this order. It partial-orders the elements of each equivalence class by attributing to each a word according to the number of its neighbours in each equivalence class: let x be a vertex; denote the number of neighbours of x in Ci by Ni(x); then x is characterised by a an ordered set of numbers (N1(x), N2(x),..., Nn(x)) which constitutes a word; we partial-order each equivalence class by the lexicographic order of the words corresponding to the vertices in the class.

Let U(X) denote the set of those vertices of X whose equivalence class is not a singleton.

If |U(X)| is bigger than a certain threshold value we jump to procedure C, which we expect to be more efficient at handling those graphs.

Procedure B: Brute force permutation procedure.

The input is a graph X together with a Canonical Vertex Classification P.

P has at most |U(X)|! refinements to a linear order.

To each of these linear orders L, compute the 0-1 array A(L) = (aij) where i belongs to U(X), j belongs to V, aij = 1 or 0 according to weather i and j are adjacent or not, and the rows and columns of the array are arranged in the order defined by L.

Let Lo be any of those L's for which A(Lo) is smaller or equal to A(L) for every L. Output Lo, end.

Any of two such Lo's differ in an automorphism of X, fixing all vertices outside U(X). So we say that Lo is a Canonical Label.

Procedure A is a breadth-first search. Procedure B is a brute force which is efficient if preceded by a very successful breadth-first search.

The following procedures constitute a depth-first search.

Procedure C:

The input is a graph X and a positive integer s.

N is the number of ways to choose s vertices out of the total number of vertices of the graph multiplied by the factorial of those s vertices (number of permutations).

This procedure generates all s-sequences of vertices of X and then applies procedure D to each of them.

If the pair (X, Si) is rejected for every i = 1,..., N, reject the pair (X, s). Else let Ai be the adjacency matrix corresponding to the relative Label with respect to Si, and find the lexicographically first among these matrices.

Output the corresponding Canonical Label.

Procedure D:

The input is a pair (X, S) where S is a sequence of vertices of X.

X/S is the complement of S in X.

To each vertex v belonging to X/S assign a binary word a(v) where aiv =1 or 0 according to weather xi and v are adjacent or not.

If these words are not all different, reject the pair (X, S).

Else order V by xi < xj if i < j; xi < v for all v belonging to X/S; and v < w for v, w belonging to X/S if a(v) < a(w).

Output (X, <).

The Canonical Label Algorithm is then:

1 - Apply A. If A rejects X, go to 3. Else

2 - Apply B to the output of A. Output the result, end.

3 - Set s = 1.

4 - Apply D to the pair (X, s). If accepted, output the result, end. Else

5 - Set s = s + 1 and go to 4.

We store all the Canonical Labels attributed to vertices contained in a basket and compare them, taking care to ensure that the vertex for which we extracted the sub-graph X1 is mapped to the corresponding vertex in the sub-graph X2.

Successive colour-refinement steps splitting vertices into classes
Figure 2 — Procedure A, step by step. G(N = 60, p = 0.05). The classification starts from valence alone and refines by neighbour-class words. Vertices sharing a colour are still indistinguishable at that step. Most classes collapse to singletons within two or three iterations — which is why the breadth-first procedure dominates the running time for almost all random graphs.

Variety Algorithm

This is the entire algorithm:

1 - Create random graph.

2 - Set Variety = 0, n = 1, basket = V (entire graph).

3 - For each element of the basket extract the number of vertices of the n_th order view (n_valence of that vertex). If n >1 & for all the elements in the basket this n_valence is the same as the (n-1)_valence then go to the next basket.

4 - If n > 1 & last basket then return to next basket in 10 (n = n -1).

5 - Subdivide the basket into smaller baskets containing vertices with the same n_valence and order the resulting baskets. Go to first basket.

6 - If the basket has only a single vertex then Variety = Variety + 1/n and go to the next basket.

7 - If n > 1 & last basket then return to next basket in 10 (n = n - 1) else if n =1 & last basket end.

8 - Extract n_th order view of each vertex belonging to the basket and apply the Canonical Label Algorithm to the sub-graph.

9 - Separate vertices with the same Canonical Label into smaller baskets and order the resulting baskets. Go to first basket.

10 - If the basket has only a single vertex then Variety = Variety + 1/n and go to the next basket.

11 - If n>1 & last basket then go to next basket in 6 else if n =1 & last basket end.

12 - Set n = n +1 and go to 3.

Considerations

The algorithm was implemented in C and it required the standard math and input/output libraries.

The decision to apply either procedure B or procedure C is made according to the fraction of |U(X)| which belongs to a single class and the size of this class.

The relevant time consuming part of the algorithm resides in the calculation of the permutations of the vertex labels.

Certain strategies were followed in the implementation in order to make the algorithm more efficient. The most relevant of those is the following:

To avoid calculating all the permutations we may consider only those that decrease the position in the lexicographic order. Since the calculation of the permutations is itself ordered, this is easily accomplished.

Analysis

For most part of the probability range the Algorithm is very efficient, since the running time is dominated by Procedure A, and so we are able to compute, in a couple of minutes, graphs with up to 1000 vertices.

For a short range of probabilities close to zero (or on the symmetric position in relation to 0.5) it might take a much longer time, due, probably, to an increase in symmetry and disconnectivity of the graph.

This range's width depends inversely on the number of vertices of the graph.

Variety plotted against edge probability p
Figure 3 — Variety against edge probability. G(N = 1000, p). Left: the full curve, symmetric about p = 0.5. Right: the transition region on a logarithmic axis. Variety rises sharply from zero and saturates near 1 across most of the range; the narrow band close to p = 0 (and its mirror near p = 1) is where symmetry and disconnectivity make the graph expensive to resolve.

References

  1. I. Ponomarenko, "Graph Algebras and the Graph Isomorphism Problem", AAECC 5, 277-286 (1994). https://doi.org/10.1007/BF01225642
  2. D. Corneil & R. Read, "The Graph Isomorphism Disease", Journal of Graph Theory, Vol. 1 (1977) 339-363. https://onlinelibrary.wiley.com/doi/abs/10.1002/jgt.3190010410
  3. L. Smolin, "Variety, complexity and cosmology", Preprint. https://arxiv.org/abs/hep-th/9203041
  4. L. Smolin, "Cosmology as a problem in critical phenomena", gr-qc/9505022 in "Complex systems and binary networks eds. R. Lopez-Pena, R. Capovilla, R. Garcia-Pelayo, H. Waelbroeck and F. Zertuche (Springer-Verlag, Berlin, 1995). https://arxiv.org/abs/gr-qc/9505022
  5. L. Smolin, "Space and Time in the Quantum Universe" in "Conceptual Problems of Quantum Gravity" ed. by A. Ashtekar and J. Stachel, (Birkhauser, Boston, 1991). https://philpapers.org/rec/SMOSAT
  6. L. Babai & L. Kucera, "Canonical Labelling of Graphs in linear average time," IEEE 2/79, 39-46. https://ieeexplore.ieee.org/document/4567999

Program

The implementation is a single C89 translation unit. It is reproduced below as it was written, with only three changes needed to make it compile under a modern toolchain: a missing #include <stdlib.h>, an explicit int return type on main(), and a corrected pointer dereference in the canonical-label routine.

Show the full C source (variety.c, ~1000 lines)
/* VARIETY */

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#define NMAX 1000

void make_graph();
void type_graph();
void quick_sort(int *b[], int left, int right);
void swap(int *b[], int i, int j);
int sub_mat(int vert);
int sub_mat_rep(int i, int j, int flag3);
void add_val(int vert);
void k_val(int ite, int lef, int rig);
void first_val();

void can_label(int ver,int iter, int pos);
void label(int nodes);
void lexico(int lef, int rig, int ki);
int divis(int i, int j, float mod);
void permut(int ver, int i, int j);
void lexico_can(int lef, int rig);
void sub_fval(int ver);

void procedure_d(int ver);
void s_sequence(int k, int level, int ver);
void permutation(int i, int j);
void procedure_c();
void same(int lef, int rig, int ki);

int a[NMAX][NMAX], val[NMAX], fval[NMAX], firval[NMAX], smat[NMAX][NMAX],
    fsmat[NMAX][NMAX], firsmat[NMAX][NMAX];
int *b[NMAX], *aux[NMAX], lote[NMAX][3];

float p;
/* double p; */
double conv, variety;
int i, j, k;
int iseed, indifference; /* iseed=random seed */

char xs;

int largo, largo2, largomajor;

int *c[NMAX], ucol[NMAX], ulin[NMAX], *d[NMAX], *e[NMAX], *e2[NMAX], **f[NMAX],
    *swp, *swdp, order[NMAX], ord[NMAX], sa[NMAX], sb[NMAX], si[NMAX], cs[NMAX],
    *h[NMAX], nodes, u, v, class, canclass, s, o, jup, vertices, klass;

int *labp, *classp, *canclassp, *permutp, *labelp, *wordp, *canlabelp;

int clip[NMAX], *clp[NMAX];

int count, first;
int inic;
int outofhere;

/*------------------------------------------------------------------*/

int main()
{
  int iter, lef, rig, one, r, runs;

  FILE *file, *fopen();

  file=fopen("var.out","a");

  iseed=55577;
  one=0;

  printf("\n\nNumber of vertices = %d\n\n", NMAX);
  fprintf(file, "\nNumber of vertices = %d\n", NMAX);

  printf("\nThe data will be written in a file called <var.out>\npress return to continue");

  xs=getchar();

  printf("\nprobability? (number between 0 and 1; if you want to exit just write 0)\n");
  scanf("%f", &p);
  printf("probability = %f\n", p);

  while(p!=0){

    printf("\n number of trials? ");
    scanf("%d", &runs);

    fprintf(file, "\nProbability = %f , number of trials = %d\n\nVariety:\n\n", p, runs);

    for(r=0; r<runs; r++){

      largomajor=0;

      indifference=0;
      variety=.0;
      conv=.0;
      lef=0;
      rig=NMAX-1;

      printf("\ninitializing variables...\n");

      for(i=0;i<NMAX;i++)
          a[i][i]=0;

      for(i=0;i<NMAX;i++)
          b[i]=&val[i];

      for(i=0; i<NMAX; i++){
          for(j=1; j<NMAX; j++)
            smat[i][j]=-1;
          lote[i][2]=1;
          lote[i][0]=lote[i][1]=0;
          smat[i][0]=i;
      }

      printf("done\n");

      srand(iseed);

      printf("\ncalculating the adjacency matrix...\n");

      make_graph();

      printf("done\n");

      if(p>.5){ /* transforms adjacency matrix into anti-adjacency matrix */
        for(i=0; i<NMAX; i++)
            for(j=0; j<NMAX; j++)
              a[i][j]=1-a[i][j];

        for(i=0;i<NMAX;i++)
            a[i][i]=0;
      }

      first_val();

      printf("\n\nnumber of vertices completed:\n");

      k_val(1,lef,rig);

      variety/=NMAX;

      printf("\n\n    VARIETY =%f    indifference =%3d    largomajor =%2d", variety,indifference,largomajor);

      printf("\nfinished trial%3d\n", r+1);

      fprintf(file, "%f\n", variety);

    }

    printf("\nprobability? (number between 0 and 1; if you want to stop just write 0) ");
    scanf("%f", &p);

  }

  printf("\nTHE END\n\n");

}
/*------------------------------------------------------------------*/

/* makes a graph with NMAX nodes and prob p */

void make_graph()
{
  int r;
  double size;

  size=32768;
  iseed=rand();
  srand(iseed);

  for(j=0;j<NMAX;j++){
    for(k=0; k<j; k++){
      r=rand();
      if((r/size) < p )
          a[j][k]=1;
      else
          a[j][k]=0;
      a[k][j]=a[j][k];
    }
  }
}

/*------------------------------------------------------------------*/

/* prints graph */

void type_graph()
{
  printf("\n N = %d , p = %f , iseed = %d \n\n  ",NMAX,p,iseed);

  for(i=0; i<NMAX; i++)
    printf("%2d", i);
  printf("\n\n");

  for(j=0;j<NMAX;j++)
    {
      printf("%2d ", j);
      for(k=0;k<NMAX;k++)
          {
            printf("%2d",a[j][k]);
          }
      printf("\n");
    }
  printf("\n");
}

/*------------------------------------------------------------------*/

/* creates matrix with first_valences */

void first_val()
{
  for(i=0; i<NMAX; i++)
    firval[i]=0;

  for(i=0; i<NMAX; i++)
    for(j=0; j<NMAX; j++)
      firsmat[i][j]=-1;

  for(i=0; i<NMAX; i++){
/*  firsmat[i][0]=i; */
    k=0;
    for(j=0; j<NMAX; j++)
      if(a[i][j]==1){
          firval[i]+=1;
          firsmat[i][k]=j;
          k++;
      }
  }

}

/*------------------------------------------------------------------*/

/* sorts val[left]...val[right] into increasing order */

void quick_sort(int *b[], int left, int right)
{
  int ii, last;

  if (left >= right)
    return;
  swap(b, left, (left+right)/2);
  last=left;
  for (ii=left+1; ii<=right; ii++)
    if (*b[ii]<*b[left])
      swap(b, ++last, ii);
  swap(b, left, last);
  quick_sort(b, left, last-1);
  quick_sort(b, last+1, right);
}

/*------------------------------------------------------------------*/

/* swap: interchange v[i] and v[j] */
void swap(int *b[], int ii, int jj)
{
  int *temp;

  temp=b[ii];
  b[ii]=b[jj];
  b[jj]=temp;
}

/*------------------------------------------------------------------*/

/* extracts the sub-matrix with the k_valences */

int sub_mat(int vert)
{
  int flg1,flg2,flag1,flag2,flag3,jmp;

  vert=b[vert]-val;

  jmp=0;

  flg2=lote[vert][2];
  flg1=lote[vert][1];
  flag1=lote[vert][0];
  k=flag2=flag3=flg2;

  for(i=flag1; i<flg2; i++)
    *aux[i]=smat[vert][i];

  while(flg1<flg2){
    i=smat[vert][flg1];
    for(j=flag2; j<flag3; j++)
      *aux[j]=smat[vert][j];

    quick_sort(aux, flag1, flag3-1);

    inic=flag1;
    for(j=0; firsmat[i][j]!=-1; j++)
      if(sub_mat_rep(i,j,flag3)==0){
          smat[vert][k]=firsmat[i][j];
          k++;
      }
    flg1++;
    flag2=flag3;
    flag3=k;

  }

  lote[vert][0]=lote[vert][1];
  lote[vert][1]=lote[vert][2];
  lote[vert][2]=k;

  if(k==flg1)
    jmp=1;

  return jmp;
}

/*------------------------------------------------------------------*/

/* avoids repeating valences */

int sub_mat_rep(int i, int j, int flag3)
{
  int rep, ii;

  if(firsmat[i][j]>*aux[flag3-1])
    return 0;

  else if(firsmat[i][j]==*aux[flag3-1])
    return 1;

  else{
    for(ii=inic; ii<flag3; ii++)
      if(firsmat[i][j]<=*aux[ii]){
          inic=ii;
          break;
      }

    if(firsmat[i][j]==*aux[ii])
      rep=1;
    else rep=0;

    return rep;
  }
}

/*------------------------------------------------------------------*/

/* sums the nodes from the k_valences */

void add_val(int vert)
{
  vert=b[vert]-val;

  for(i=1; smat[vert][i]!=-1 && i<NMAX; i++);

  val[vert]=i-1;
}

/*------------------------------------------------------------------*/

/* recursively calculates the k_valences, orders the vertices and calculates the variety */

void k_val(int ite, int lef, int rig)
{
  int flg1, flg2, jump, vertex, fg1, fg2, ver;

  jump=0;

  printf("\n(calculating sub_matrix; iteration # %d\n",ite);
  for(vertex=lef; vertex<=rig; vertex++){

    for(i=0; i<NMAX; i++)
      aux[i]=&smat[b[vertex]-val][i];

    jump=sub_mat(vertex);
    add_val(vertex);
  }
  printf("         ... done)");

  quick_sort(b,lef,rig);

  flg1=flg2=lef;

  while(flg2<=rig){

    printf("\n %d",flg2);

    while(flg2<rig && *b[flg2]==*b[flg2+1])
      flg2++;

    if(flg1==flg2){
      indifference+=ite;
      conv=ite;
      variety+=1.0/conv;
      jump=0;
      flg2++;
      flg1=flg2;
    }

    else if(jump==1){  /* eliminates the zero valences */

      jump=0;
      flg2++;

      flg1=flg2;
    }

    else{

      nodes=val[b[flg1]-val]+1;
      vertices=flg2-flg1+1;

      canlabelp=(int *)calloc(nodes*nodes*vertices, sizeof(int));

      for(ver=flg1; ver<=flg2; ver++){

          labelp=(int *)calloc(nodes, sizeof(int));

          classp=(int *)calloc(nodes, sizeof(int));

          can_label(b[ver]-val,ite, ver); /* stores a canonical label in labelp */

/* convert canlabel into number -- DEAD CODE, commented out
p=0;
q=r=1;
  for(j=nodes-1; j>=0; j--)
    for(k=nodes-1; k>=0; k--){
      if(q>nodes*nodes) goto out;
      else if(r<128) canlabel+=pow(2,r)*a[*(labelp+j)][*(labelp+k)];

      else{
            r=0;
            *(canlabelp+(ver-flg1)*(nodes*nodes)/127+p)=canlabel;
            canlabel=0;
            p++;
      }
      q++;
      r++;
    }
  out:;
*/

          for(j=0; j<nodes; j++)
            for(k=0; k<nodes; k++)
              *(canlabelp+(ver-flg1)*nodes*nodes+j*nodes+k)=a[*(labelp+j)][*(labelp+k)];


          free(classp);

          free(labelp);

      }

      canclass=0;
      k=0;

      canclassp=(int *)calloc(vertices, sizeof(int));

      lexico_can(0,vertices-1);

      free(canlabelp);

      for(i=flg1; i<=flg2; i++)
          *b[i]=*(canclassp+i-flg1);

      free(canclassp);

      quick_sort(b, flg1, flg2);

      fg1=fg2=flg1;
      while(fg2<=flg2){

          while(fg2<flg2 && *b[fg2]==*b[fg2+1])
            fg2++;

          if(fg1==fg2){
            indifference+=ite;
            conv=ite;
            variety+=1.0/conv;
            fg2+=1;
            fg1=fg2;
          }

          else{
            k_val(ite+1,fg1,fg2);
            fg2+=1;
            fg1=fg2;
          }
      }
    }
  flg2+=1;
  flg1=flg2;
  }

return;

}

/*------------------------------------------------------------------*/
/*------------------------------------------------------------------*/

/* atributes a canonical label to a graph */

void can_label(int ver,int iter, int pos)
{
  int lg;

  sub_fval(ver);

  labp=(int *)calloc(nodes*(nodes+1), sizeof(int));
  label(nodes);
  free(labp);

  if(u<2){

    for(i=0; i<nodes; i++)
      *(labelp+i)=smat[ver][**f[i]];
  }

  else{

    lg=log(nodes)/log(log(nodes));

    if(largo>9 && largo>u*2.0/3.0){       /* (u > lg)  7! = 5040 , 8! = 40320 */

      printf("\n working on %d th vertex , largo = %d , u = %d", pos, largo, u);
      printf("\n TROUBLE: (10! = 3628800 trials) If our breadth-first procedure is not efficient enough there is a real possibility that this might take quite a long while because we might have to calculate all %d! permutations\n", u);

      procedure_d(ver);

    }

    else{

      if(largo > 8 || (largo>7 && largo2>1))
          printf("\n\n This might take a while since we have to calculate something of the order of %d! * %d! (9! = 362880) , u = %d\n", largo, largo2, u);

      permutp=(int *)calloc(u*nodes, sizeof(int));
      swp=(int *)calloc((u-1)*2, sizeof(int));

      for(i=0; i<u; i++)
          for(j=0; j<nodes; j++)
            *(permutp+i*nodes+j)=a[smat[ver][*d[i]]][smat[ver][**f[j]]];

      for(i=0; i<nodes; i++)
          *(labelp+i)=smat[ver][**f[i]];

      permut(ver,0,0);

      free(swp);
      free(permutp);

    }
  }
}

/*------------------------------------------------------------------*/

/* labelling a graph */

void label(int nodes)
{
  float mod;

  mod=4.0;

  for(i=0; i<nodes; i++)
    for(j=0; j<nodes; j++)
      *(labp+i*nodes+j)=0;

  for(i=0; i<nodes; i++){
    clp[i]=(labp+i*nodes);
    *(labp+i*nodes)=fval[i];
  }

  quick_sort(clp,0,nodes-1);

  klass=1;
  clip[0]=klass;
  for(i=0; i<nodes-1; i++){
    if(*clp[i]!=*clp[i+1]) klass++;
    clip[i+1]=klass;
  }

  for(i=0; i<nodes; i++)
    *clp[i]=clip[i];

  for(i=0; i<nodes; i++)
    for(j=0; fsmat[i][j]!=-1; j++)
      (*(labp+i*nodes+*(labp+fsmat[i][j]*nodes)))++;

/* we could do mod4 */
/* for(i=0; i<nodes; i++)
    for(j=1; j<=klass; j++)
      *(labp+i*nodes+j)=divis(i,j,mod);
*/

/* convert label into number -- DEAD CODE, commented out
  p=q=0;

  while(q<nodes){
    for(i=128*2/mod-1; i>=0; i--){
        if(q<nodes) label+=pow(mod,i)*(*(labip+q));
        else break;
        q++;
    }
    *(labp+p)=label;
    p++;
  }
*/

  largo2=largo=0;
  u=v=0;
  class=0;

  lexico(0, nodes-1, 0);

}

/*------------------------------------------------------------------*/

/* lexicographic order */

void lexico(int lef, int rig, int ki)
{
  int flg1,flg2; /* *c[nodes], ucol[nodes], ulin[nodes] */

  if(ki==0)
    for(i=0; i<nodes; i++){
      c[i]=(labp+i*nodes);
      ulin[i]=-1;
      ucol[i]=-1;
      *(classp+i)=-1;
    }

  else
    for(i=lef; i<=rig; i++)
      c[i]++;

  quick_sort(c,lef,rig);

  flg1=flg2=lef;

  while(flg2<=rig){
    while(flg2<rig && *c[flg2]==*c[flg2+1])
      flg2++;

  if(flg1!=flg2 && ki<klass) lexico(flg1,flg2, ki+1);

  else{
    for(i=flg1; i<=flg2; i++){

        if(flg2-flg1+1>largo){
          largo2=largo;
          largo=(flg2-flg1+1);
        }

        if(flg2-flg1+1>largomajor)
          largomajor=(flg2-flg1+1);

        ucol[i]=(c[i]-labp)/nodes;
        *(classp+ucol[i])=class;
        if(flg1!=flg2){
          ulin[u]=ucol[i];
          d[u]=&ulin[u];
          f[i]=&e[u];
          e[u]=&ucol[i];
          u++;
        }
        else{
          f[i]=&e2[v];
          e2[v]=&ucol[i];
          v++;
        }
    }
    class++;
  }
  flg2+=1;
  flg1=flg2;
  }
}

/*------------------------------------------------------------------*/

/* mod function */

int divis(int i, int j, float mod)
{
  int x;

  x=*(labp+i*nodes+j)/mod;
  x=(*(labp+i*nodes+j)/mod-x)*mod;

  return x;
}

/*------------------------------------------------------------------*/

/* permutations: returns canonical label in *f[nodes] */

void permut(int ver, int i, int j)
{
  int equal;

  while(j<u-1){
    while(j<u && *(classp+*d[i])==*(classp+*d[j])){

      *(swp+i)=i;
      *(swp+u-1+i)=j;

      permut(ver,i+1,i+1);
      j++;
    }
    return;
  }

  for(i=0; i<u-1; i++){
    swap(d,*(swp+i),*(swp+u-1+i));
    swap(e,*(swp+i),*(swp+u-1+i));
  }


  equal=1;
  for(i=0; i<u; i++){
    for(j=0; j<nodes; j++){
      if(equal==1){
          if(a[smat[ver][*d[i]]][smat[ver][**f[j]]]>*(permutp+i*nodes+j)) goto out;
          else if(a[smat[ver][*d[i]]][smat[ver][**f[j]]]<*(permutp+i*nodes+j)) equal=0;
      }
      *(permutp+i*nodes+j)=a[smat[ver][*d[i]]][smat[ver][**f[j]]];
    }
  }
  for(i=0; i<nodes; i++)
    *(labelp+i)=smat[ver][**f[i]];

  out:;

  for(i=u-2; i>=0; i--){
    swap(d,*(swp+i),*(swp+u-1+i));
    swap(e,*(swp+i),*(swp+u-1+i));
  }

  return;

}

/*------------------------------------------------------------------*/

/* lexicographic order of canonical labels */

void lexico_can(int lef, int rig)
{
  int flg1,flg2;

  if(k==0)
    for(i=0; i<vertices; i++)
      c[i]=(canlabelp+i*nodes*nodes);

  else
    for(i=lef; i<=rig; i++)
      c[i]+=1;

  k++;

  quick_sort(c,lef,rig);

  flg1=flg2=lef;

  while(flg2<=rig){
    while(flg2<rig && *c[flg2]==*c[flg2+1])
      flg2++;

    if(flg1!=flg2 && k<nodes*nodes) lexico_can(flg1,flg2);

    else{
      for(i=flg1; i<=flg2; i++)
          *(canclassp+(c[i]-canlabelp)/(nodes*nodes))=canclass;
      canclass++;
    }

    flg2+=1;
    flg1=flg2;
  }

  return;

}

/*------------------------------------------------------------------*/

/* creates matrix with first_valences of sub_graph */

void sub_fval(int ver)
{
  for(i=0; i<nodes; i++)
    fval[i]=0;

  for(i=0; i<nodes; i++)
    for(j=0; j<nodes; j++)
      fsmat[i][j]=-1;

  for(i=0; i<nodes; i++){
    k=0;
    for(j=0; fsmat[ver][j]!=-1 && j<nodes; j++)
      if(a[smat[ver][i]][smat[ver][j]]==1){
          fval[i]+=1;
          fsmat[i][k]=j;
          k++;
      }
  }

}

/*------------------------------------------------------------------*/
/*------------------------------------------------------------------*/

/* tests the difficult graphs */
/* outputs canonical label in *labelp[nodes] */

void procedure_d(int ver)
{
/* s=6*sqrt(nodes); */
/* s=sqrt(nodes); */

  s=1;
  o=0;

  while(s<nodes){

    count=0;

    swdp=(int *)calloc(s, sizeof(int));
    s_sequence(0,0,ver);
    free(swdp);

    if(o!=0)
      break;

    else
      s++;
  }

}

/*------------------------------------------------------------------*/

/* extracts s elements from smat[nodes] */

void s_sequence(int k, int level, int ver)
{
  int ii, kk;

  while(level<s){

    for(ii=k; ii<=(level+nodes-s); ii++){
      sa[level]=smat[ver][ii];
      si[level+1]=ii;
      s_sequence(ii+1, level+1, ver);
    }
    return;
  }
  si[0]=-1;
  si[level+1]=nodes;

  kk=0;
  for(i=0; i<s+1; i++) /* vertices not in the s-sequence */
    if(si[i]+1!=si[i+1])
      for(j=si[i]+1; j<si[i+1]; j++){
          cs[kk]=smat[ver][j];
          kk++;
      }

  first=1;

  permutation(0,0);

  return;

}

/*------------------------------------------------------------------*/

/* permutes the s elements from smat[nodes] */

void permutation(int i, int j)
{
  int temp, equal, ii, jj;

  jup=0;

  while(j<s-1){
    while(j<s){

      *(swdp+i)=i;
      *(swdp+s-1+i)=j;

      permutation(i+1, i+1);
      j++;

    }
    return;
  }

  for(ii=0; ii<s; ii++)
    sb[ii]=sa[ii];

  for(ii=0; ii<s-1; ii++){
    temp=sb[*(swdp+ii)]; /* swap(i,j) */
    sb[*(swdp+ii)]=sb[*(swdp+s-1+ii)];
    sb[*(swdp+s-1+ii)]=temp;
  }

  count++;

  wordp=(int *)calloc((nodes-s)*s, sizeof(int));
  procedure_c();
  free(wordp);

  if(jup==0){

    if(first==1){
      for(i=0; i<s; i++)
          *(labelp+i)=sb[i];
      for(i=s; i<nodes; i++)
          *(labelp+i)=ord[i-s];
      first=0;
    }


    equal=1;
    for(j=0; j<nodes; j++)
      for(k=0; k<nodes; k++) /* canonical label in *labelp */
          if(k<s){
            if(equal==1){
              if(a[sb[j]][sb[k]]>a[*(labelp+j)][*(labelp+k)]) goto out;
              else if (a[sb[j]][sb[k]]<a[*(labelp+j)][*(labelp+k)]) equal=0;
            }
            *(labelp+k)=sb[k];
          }
          else{
            if(equal==1){
              if(a[ord[j-s]][ord[k-s]]>a[*(labelp+j)][*(labelp+k)]) goto out;
              else if (a[ord[j-s]][ord[k-s]]<a[*(labelp+j)][*(labelp+k)]) equal=0;
            }
            *(labelp+k)=ord[k-s];
          }

  out:;

    o++;
  }
  else jup=0;

  return;

}

/*------------------------------------------------------------------*/

/* word[i][j] */

void procedure_c() /* wordp */
{
  for(i=0; i<nodes-s; i++){
    for(j=0; j<s; j++){
      *(wordp+i*s+j)=a[cs[i]][sb[j]];
    }
  }

  outofhere=0;

  same(0, nodes-s-1, 0);

}

/*------------------------------------------------------------------*/

/* checks if words are all different and orders them */

void same(int lef, int rig, int ki) /* ord */
{
  int flg1,flg2; /* *h[nodes-s] */

  if(outofhere==1)
    return;

  if(ki==0)
    for(i=0; i<nodes-s; i++)
      h[i]=(wordp+i*s);

  else
    for(i=lef; i<=rig; i++)
      h[i]+=1;

  quick_sort(h,lef,rig);

  flg1=flg2=lef;

  while(flg2<=rig){
    while(flg2<rig && *h[flg2]==*h[flg2+1])
      flg2++;

  if(flg1!=flg2 && ki<s-1)
    same(flg1,flg2,ki+1);

  else{
    for(i=flg1; i<=flg2; i++){

        if(flg1!=flg2){
          jup=1;
          outofhere=1;
          goto same_word; /* tenho que sair destes loops */
        }

        else{
          ord[i]=cs[(h[i]-wordp)/s];
        }
      }
    }
  flg2+=1;
  flg1=flg2;
  }

  same_word:;

}

/*------------------------------------------------------------------*/

/* THE END */

Source, figures and a modern WL-based reimplementation: github.com/ferrelm/variety-graph

Comentários

Mensagens populares deste blogue

ITRA Performance Index - Everything You Always Wanted to Know But Were Afraid to Ask

Provas Insanas - Westfield Sydney to Melbourne Ultramarathon 1983

The Ministry of Doubt

Andorra Ultra Trail VallNord - Ronda dels Cims 2019 - Prólogo

The Unreliable Agent: Why Guardrails Are Not Guarantees

5º Ultra Trail da Serra da Freita 2010 - Recordações