Schneller Algo + Berechnung des Vereinfachten Bruhat-Graphen

Makefile

@@ -1,7 +1,7 @@
 HEADERS=coxeter.h thickenings.h queue.h
 OPTIONS=-O3 -g -std=gnu99
 
-all: generate process
+all: generate process test
 
 generate: generate.o coxeter.o thickenings.o
 	gcc $(OPTIONS) -o generate generate.o thickenings.o coxeter.o -lgsl -lcblas
@@ -9,6 +9,9 @@ generate: generate.o coxeter.o thickenings.o
 process: process.o coxeter.o thickenings.o
 	gcc $(OPTIONS) -o process process.o thickenings.o coxeter.o -lgsl -lcblas
 
+test: test.o coxeter.o thickenings.o
+	gcc $(OPTIONS) -o test test.o thickenings.o coxeter.o -lgsl -lcblas
+
 generate.o: generate.c $(HEADERS)
 	gcc $(OPTIONS) -c generate.c
 
@@ -21,5 +24,8 @@ thickenings.o: thickenings.c $(HEADERS)
 coxeter.o: coxeter.c $(HEADERS)
 	gcc $(OPTIONS) -c coxeter.c
 
+test.o: test.c $(HEADERS)
+	gcc $(OPTIONS) -c test.c
+
 clean:
-	rm -f generate thickenings.o coxeter.o generate.o process.o
+	rm -f generate process test thickenings.o coxeter.o generate.o process.o test.o

process.c

@@ -47,7 +47,8 @@ int main(int argc, const char *argv[])
   int rank, order;
   semisimple_type_t type;
   int n;
-  int *thickenings, *level;
+  int *thickenings;
+  signed char *level;
   node_t *graph;
   int *left, *right;
   int *left_invariant, *right_invariant;
@@ -58,7 +59,7 @@ int main(int argc, const char *argv[])
   int current;
   int *seen;
   int *generators;
-  int ngens;
+  int ngens, nleft, nright;
 
   char string_buffer1[1000];
   const char *alphabet = "abcdefghijklmnopqrstuvwxyz";
@@ -86,7 +87,7 @@ int main(int argc, const char *argv[])
   right_invariant = (int*)malloc(rank*sizeof(int));
   seen = (int*)malloc(order*sizeof(int));
   generators = (int*)malloc(order*sizeof(int));
-  level = (int*)malloc(order*sizeof(int));
+  level = (signed char*)malloc(order*sizeof(int));
   for(int i = 0; i < order; i++) {
     graph[i].left = &left[i*rank];
     graph[i].right = &right[i*rank];
@@ -97,11 +98,13 @@ int main(int argc, const char *argv[])
 
   int counter = 0;
 
-  while(fread(level, sizeof(int), order, infile) == order) {
+  while(fread(level, sizeof(signed char), order, infile) == order) {
 
+    /*
     if((counter++) % 100000 == 0)
       print_thickening(rank, order, level, 0, 0, 0, alphabet, stdout);
     continue;
+    */
 
     for(int j = 0; j < rank; j++) {
       left_invariant[j] = 1;
@@ -109,17 +112,27 @@ int main(int argc, const char *argv[])
     }
     for(int j = 0; j < order; j++) {
       for(int k = 0; k < rank; k++) {
-	if(level[j] > 0 && level[graph[j].left[k]] < 0 || level[j] < 0 && level[graph[j].left[k]] > 0)
+	if(level[j] > 0 && level[graph[j].left[k]] < 0 || level[j] < 0 && level[graph[j].left[k]] > 0) {
 	  left_invariant[k] = 0;
-	if(level[j] > 0 && level[graph[j].right[k]] < 0 || level[j] < 0 && level[graph[j].right[k]] > 0)
+	}
+	if(level[j] > 0 && level[graph[j].right[k]] < 0 || level[j] < 0 && level[graph[j].right[k]] > 0) {
 	  right_invariant[k] = 0;
+	}
       }
     }
 
+    nleft = nright = 0;
+    for(int j = 0; j < rank; j ++) {
+      if(left_invariant[j])
+	nleft++;
+      if(right_invariant[j])
+	nright++;
+    }
+
     ngens = 0;
     memset(generators, 0, order*sizeof(int));
     for(int j = 0; j < order; j++) {
-      if(level[j] == 1000 && generators[j] == 0) { // ignore the generator, if it is equivalent to one already seen
+      if(level[j] == HEAD_MARKER && generators[j] == 0) { // ignore the generator, if it is equivalent to one already seen
 	ngens++;
 	queue_init(&queue);
 	queue_put(&queue, j);
@@ -137,31 +150,33 @@ int main(int argc, const char *argv[])
       }
     }
 
-    printf("left: ");
+    if(nleft >= 2 && nright >= 2) {
-    for(int j = 0; j < rank; j++)
+      printf("left: ");
-      printf("%c", left_invariant[j] ? alphabet[j] : ' ');
+      for(int j = 0; j < rank; j++)
-    printf(" right: ");
+	printf("%c", left_invariant[j] ? alphabet[j] : ' ');
-    for(int j = 0; j < rank; j++)
+      printf(" right: ");
-      printf("%c", right_invariant[j] ? alphabet[j] : ' ');
+      for(int j = 0; j < rank; j++)
-    printf(" generators: ");
+	printf("%c", right_invariant[j] ? alphabet[j] : ' ');
+      printf(" generators: ");
 
-    memset(seen, 0, order*sizeof(int));
+      memset(seen, 0, order*sizeof(int));
-    for(int i = 0; i < order; i++) {
+      for(int i = 0; i < order; i++) {
-      if(generators[i] != 0 && seen[generators[i]-1] == 0) {
+	if(generators[i] != 0 && seen[generators[i]-1] == 0) {
-	seen[generators[i]-1] = 1;
+	  seen[generators[i]-1] = 1;
-	//	if(type.n == 1 && type.factors[0].series == 'A')
+	  //	if(type.n == 1 && type.factors[0].series == 'A')
-	//	  printf("%s ", stringify_SLn1_permutation(graph[i].word, graph[i].wordlength, rank, string_buffer1));
+	  //	  printf("%s ", stringify_SLn1_permutation(graph[i].word, graph[i].wordlength, rank, string_buffer1));
-	//	else if(type.n == 1 && (type.factors[0].series == 'B' || type.factors[0].series == 'C'))
+	  //	else if(type.n == 1 && (type.factors[0].series == 'B' || type.factors[0].series == 'C'))
-	//	  printf("%s ", stringify_Onn1_permutation(graph[i].word, graph[i].wordlength, rank, string_buffer1));
+	  //	  printf("%s ", stringify_Onn1_permutation(graph[i].word, graph[i].wordlength, rank, string_buffer1));
-	//	else
+	  //	else
 	  if(i == 0)
 	    printf("1 ");
 	  else
 	    printf("%s ", alphabetize(graph[i].word, graph[i].wordlength, alphabet, string_buffer1));
+	}
       }
-    }
 
-    printf("\n");
+      printf("\n");
+    }
   }
 
   if(infile != stdin)

test.c

@@ -0,0 +1,97 @@
+#include <stdio.h>
+
+#include "thickenings.h"
+#include "queue.h"
+
+int main(int argc, const char *argv[])
+{
+  semisimple_type_t type;
+  node_t *graph;
+  int *leftbuf, *rightbuf;
+  edgelist_t *edgelists;
+  int *words;
+  int rank, order;
+  int *reduced, *group;
+  int current;
+  queue_t queue;
+  char alphabet[] = "abcdefghijklmnopqrstuvwxyz";
+  char buffer[1024];
+
+  int left = ~0x01; // 1111 1110
+  int right = ~0x01; // 1111 1101
+
+  type.n = 1;
+  type.factors = (simple_type_t*)malloc(type.n*sizeof(simple_type_t));
+  type.factors[0].series = 'B';
+  type.factors[0].rank = 4;
+
+  rank = coxeter_rank(type);
+  order = coxeter_order(type);
+  graph = (node_t*)malloc(order*sizeof(node_t));
+  leftbuf = (int*)malloc(rank*order*sizeof(int));
+  rightbuf = (int*)malloc(rank*order*sizeof(int));
+  for(int i = 0; i < order; i++) {
+    graph[i].left = &leftbuf[i*rank];
+    graph[i].right = &rightbuf[i*rank];
+  }
+
+  prepare_graph(type, graph, &edgelists, &words);
+
+  reduced = (int*)malloc(order*sizeof(int));
+  group = (int*)malloc(order*sizeof(int));
+  for(int i = 0; i < order; i++) {
+    group[i] = -1;
+    reduced[i] = i;
+  }
+
+  // step 1: group
+  for(int i = 0; i < order; i++) {
+    if(group[i] != -1)
+      continue;
+
+    queue_init(&queue);
+    queue_put(&queue, i);
+    while((current = queue_get(&queue)) != -1) {
+      if(group[current] != -1)
+	continue;
+      group[current] = i;
+
+      for(int j = 0; j < rank; j++) {
+	if(left & (1 << j))
+	  queue_put(&queue, graph[current].left[j]);
+	if(right & (1 << j))
+	  queue_put(&queue, graph[current].right[j]);
+      }
+    }
+  }
+
+  // step 2: find minimum
+  for(int i = 0; i < order; i++)
+    if(graph[i].wordlength < graph[reduced[group[i]]].wordlength)
+      reduced[group[i]] = i;
+
+  // step 3: assign minimum to all
+  for(int i = 0; i < order; i++)
+    reduced[i] = reduced[group[i]];
+
+  for(int i = 0; i < order; i++)
+    if(reduced[i] == i) {
+      if(i == 0)
+	printf("1 ");
+      else
+	printf("%s ", alphabetize(graph[i].word, graph[i].wordlength, alphabet, buffer));
+    }
+
+  printf("\n");
+
+  free(type.factors);
+  free(graph);
+  free(edgelists);
+  free(words);
+  free(reduced);
+  free(group);
+  free(leftbuf);
+  free(rightbuf);
+
+  return 0;
+}

thickenings.c

@@ -20,15 +20,15 @@ char *alphabetize(int *word, int len, const char *alphabet, char *buffer)
   return buffer;
 }
 
-void print_thickening(int rank, int order, const int *thickening, int is_fat, int is_slim, int conflict, const char *alphabet, FILE *f)
+void print_thickening(int rank, int order, const signed char *thickening, int is_fat, int is_slim, int conflict, const char *alphabet, FILE *f)
 {
   for(int i = 0; i < order; i++) {
-    if(thickening[i] == 1000)
+    if(thickening[i] == HEAD_MARKER)
-      fprintf(f, "\e[41mx\e[40m\e[m");
+      fprintf(f, "\e[41mx\e[39m\e[49m");
     else if(thickening[i] < 0 && thickening[i] > -10)
-      fprintf(f, "\e[47m\e[30m%d\e[40m\e[m", -thickening[i]);
+      fprintf(f, "\e[47m\e[30m%d\e[39m\e[49m", -thickening[i]);
     else if(thickening[i] <= -10)
-      fprintf(f, "\e[47m\e[30m+\e[40m\e[m");
+      fprintf(f, "\e[47m\e[30m+\e[39m\e[49m");
     else if(thickening[i] > 0 && thickening[i] < 10)
       fprintf(f, "%d", thickening[i]);
     else if(thickening[i] >= 10)
@@ -44,10 +44,10 @@ void print_thickening(int rank, int order, const int *thickening, int is_fat, in
   if(conflict)
     fprintf(f, " C");
 
-  fprintf(f, "\n");
+  fprintf(f, "\e[K\n");
 }
 
-void print_balanced_thickening(int rank, int order, const int *thickening, const int *left_invariant, const int *right_invariant, const char *alphabet, FILE *f)
+void print_balanced_thickening(int rank, int order, const signed char *thickening, const int *left_invariant, const int *right_invariant, const char *alphabet, FILE *f)
 {
   for(int i = 0; i < order; i++) {
     if(thickening[i])
@@ -309,13 +309,11 @@ void prepare_graph(semisimple_type_t type, node_t *graph, edgelist_t **edgelists
 void enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, const char *alphabet, FILE *outfile)
 {
   int rank, order;
-  int *level;
+  signed char *level;
-  int *left_invariant, *right_invariant;
+  long count;
-  long thickenings_count, fat_count, slim_count, balanced_count;
   int is_fat, is_slim;
   int current_level, head, current;
   int i;
-  int conflict;
   edgelist_t *edge;
 
   queue_t queue;
@@ -323,120 +321,79 @@ void enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, const
   rank = coxeter_rank(type);
   order = coxeter_order(type);
 
-  level = (int*)malloc(order*sizeof(int));
+  level = (signed char*)malloc(order*sizeof(int));
-  left_invariant = (int*)malloc(rank*sizeof(int));
-  right_invariant = (int*)malloc(rank*sizeof(int));
 
-  thickenings_count = fat_count = slim_count = balanced_count = 0;
+  count = 0;
   memset(level, 0, order*sizeof(int));
   current_level = 1;
   head = order - 1;
-  level[head] = 1000;
+  level[head] = HEAD_MARKER;
   while(current_level > 0) {
-    // calculate transitive closure
+
-    conflict = 0;
+    // calculate transitive closure; that is, fill current_level in every spot which must be marked with the current hat (but was not already before), and -current_level in every opposite spot (including opposite to the head)
+    is_slim = 1;
     queue_init(&queue);
+    level[graph[head].opposite] = -current_level;
     queue_put(&queue, head);
-    for(int i = head + 1; level[i] != 1000 && i < order; i++) {
+    for(int i = head + 1; level[i] != HEAD_MARKER && i < order; i++) { // everything which is right to the head and empty will not get marked in this or higher levels, so we can mark its opposite
-      if(level[graph[i].opposite] == 0) {
+      if(level[i] == current_level) {
+	is_slim = 0;
+	break;
+      } if(level[i] == 0) {
+	level[i] = -current_level;
 	level[graph[i].opposite] = current_level;
 	queue_put(&queue, graph[i].opposite);
       }
     }
-    while((current = queue_get(&queue)) != -1) {
+
-      if(level[current] < 0 || level[graph[current].opposite] > 0) { // conflict, can not be slim
+    if(is_slim) {
-	conflict = 1;
+      while((current = queue_get(&queue)) != -1) {
-	break;
+	edge = graph[current].bruhat_lower;
+	while(edge) {
+	  if(level[edge->to] < 0) {
+	    is_slim = 0;
+	    break;
+	  } else if(level[edge->to] == 0) {
+	    level[edge->to] = current_level;
+	    level[graph[edge->to].opposite] = -current_level;
+	    queue_put(&queue, edge->to);
+	  }
+	  edge = edge->next;
+	}
       }
-      if(level[graph[current].opposite] == 0)
+    }
-	 level[graph[current].opposite] = - current_level;
+
-      edge = graph[current].bruhat_lower;
+    // now we have something, check if it is a balanced thickening
-      while(edge) {
+    if(is_slim) {
-	if(level[edge->to] < 0) {
+      is_fat = 1;
-	  conflict = 1;
+      for(int i = 0; i < order; i++) {
+	if(level[i] == 0) {
+	  is_fat = 0;
 	  break;
 	}
-	if(level[edge->to] == 0) {
+      }
-	  level[edge->to] = current_level;
+
-	  queue_put(&queue, edge->to);
+      if(is_fat) {
-	}
+	//      ERROR(count >= MAX_THICKENINGS, "Too many balanced thickenings! Increase MAX_THICKENINGS\n");
-	edge = edge->next;
+	if(count % 10000000 == 0)
+	  print_thickening(rank, order, level, 0, 0, 0, alphabet, stderr);
+	//	fwrite(level, sizeof(signed char), order, outfile);
+	count++;
       }
     }
 
-    /*
-    // we have a thickening, do something with it!
-    is_fat = is_slim = 1;
-    for(int i = 0; i < order; i++) {
-      if(level[graph[i].opposite] != 0) {
-	if(level[i] != 0)
-	  is_slim = 0;
-      } else {
-	if(level[i] == 0)
-	  is_fat = 0;
-      }
-    }
-
-    // count
-    thickenings_count++;
-    if(is_fat)
-      fat_count++;
-    if(is_slim)
-      slim_count++;
-    if(is_slim && is_fat) {
-      ERROR(balanced_count >= MAX_THICKENINGS, "Too many balanced thickenings! Increase MAX_THICKENINGS\n");
-      //memcpy(&balanced_thickenings[balanced_count*order], level, order*sizeof(int));
-      fwrite(level, sizeof(int), order, outfile);
-      balanced_count++;
-    }
-    */
-
-    /*
-    // print out the thickening
-    if(is_fat && is_slim) {
-      // check for invariances
-      for(int j = 0; j < rank; j++) {
-    	left_invariant[j] = 1;
-    	right_invariant[j] = 1;
-      }
-      for(int i = 0; i < order; i++) {
-    	for(int j = 0; j < rank; j++) {
-    	  if(level[i] == 0 && level[graph[i].left[j]] != 0 || level[i] != 0 && level[graph[i].left[j]] == 0)
-    	    left_invariant[j] = 0;
-    	  if(level[i] == 0 && level[graph[i].right[j]] != 0 || level[i] != 0 && level[graph[i].right[j]] == 0)
-    	    right_invariant[j] = 0;
-    	}
-      }
-      print_balanced_thickening(rank, order, level, left_invariant, right_invariant, alphabet, stderr);
-    }
-    */
-
-    is_fat = 1;
-    for(int i = 0; i < order; i++) {
-      if(level[i] == 0)
-	is_fat = 0;
-    }
-
-    if(is_fat && !conflict) {
-      //      ERROR(balanced_count >= MAX_THICKENINGS, "Too many balanced thickenings! Increase MAX_THICKENINGS\n");
-      fwrite(level, sizeof(int), order, outfile);
-      balanced_count++;
-    }
-
-    // print_thickening(rank, order, level, is_fat, !conflict, conflict, alphabet, stderr);
-
     // now find the next one!
 
     // try to find empty spot to the left of "head", but only if it is slim, as otherwise there is no point in adding even more
-    if(!conflict) {
+    if(is_slim) {
       for(i = head - 1; i >= 0; i--)
 	if(level[i] == 0)
 	  break;
       if(i >= 0) {
 	head = i;
-	level[head] = 1000;
+	level[head] = HEAD_MARKER;
 	current_level++;
-	// print_thickening(rank, order, level, 0, 0, 0, alphabet, stderr);
+	ERROR(current_level >= HEAD_MARKER, "HEAD_MARKER to small!\n");
 	continue;
       }
     }
@@ -449,31 +406,25 @@ void enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, const
       if(i >= 0) { // if this was successful, just move head
 	level[head] = 0;
 	head = i;
-	level[head] = 1000;
+	level[head] = HEAD_MARKER;
 	break;
       } else { // if moving the head is not possible, take the next head to the right
 	current_level--;
 	level[head] = 0;
 	do {
 	  head++;
-	} while(head < order && level[head] != 1000);
+	} while(head < order && level[head] != HEAD_MARKER);
       }
     }
 
     // clean up
     for(int i = 0; i < order; i++)
-      if(level[i] >= current_level && level[i] != 1000 || level[i] <= -current_level)
+      if(level[i] >= current_level && level[i] != HEAD_MARKER || level[i] <= -current_level)
 	level[i] = 0;
-
-
-    //    print_thickening(rank, order, level, 0, 0, 0, alphabet, stderr);
-
   }
 
   fprintf(stderr, "\n");
-  fprintf(stderr, "Found %ld thickenings, %ld fat, %ld slim, %ld balanced\n\n", thickenings_count, fat_count, slim_count, balanced_count);
+  fprintf(stderr, "Found %ld balanced thickenings\n\n", count);
 
   free(level);
-  free(left_invariant);
-  free(right_invariant);
 }

thickenings.h

@@ -6,6 +6,7 @@
 #define DEBUG(msg, ...) do{fprintf(stderr, msg, ##__VA_ARGS__); }while(0)
 
 #define MAX_THICKENINGS 10000000
+#define HEAD_MARKER 127
 
 typedef struct _edgelist {
   int to;
@@ -24,8 +25,8 @@ typedef struct {
 } node_t;
 
 char *alphabetize(int *word, int len, const char *alphabet, char *buffer);
-void print_balanced_thickening(int rank, int order, const int *thickening, const int *left_invariant, const int *right_invariant, const char *alphabet, FILE *f);
+void print_balanced_thickening(int rank, int order, const signed char *thickening, const int *left_invariant, const int *right_invariant, const char *alphabet, FILE *f);
-void print_thickening(int rank, int order, const int *thickening, int fat, int slim, int conflict, const char *alphabet, FILE *f);
+void print_thickening(int rank, int order, const signed char *thickening, int fat, int slim, int conflict, const char *alphabet, FILE *f);
 static int compare_wordlength(const void *a, const void *b, void *gr);
 void prepare_graph(semisimple_type_t type, node_t *graph, edgelist_t **edgelists_pointer, int **words_pointer);
 void enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, const char *alphabet, FILE *outfile);