florian/enumerate-balanced-ideals
1
0
Fork 0
Code Issues Pull Requests Releases Activity

Getting better, but still not there

Browse Source
This commit is contained in:
Florian Stecker 2016-11-16 23:46:37 +01:00
parent 1811cf143e
commit ed63dc2b82
2 changed files with 30 additions and 240 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

46
bitvec.h
View File

@ -114,18 +114,17 @@ static inline int bv_get_bit(const bitvec_t *x, int k)
static inline void bv_clear(bitvec_t *x) static inline void bv_clear(bitvec_t *x)
{ {
int i; int i;
for (i=0;i<BV_QWORD_RANK;i++) { for (i=0;i<BV_QWORD_RANK;i++)
x->v[i] = 0; x->v[i] = 0;
}
} }
static inline int bv_is_zero(const bitvec_t *x) static inline int bv_is_zero(const bitvec_t *x)
{ {
int i; int i;
for (i=0;i<BV_QWORD_RANK;i++) { for (i=0;i<BV_QWORD_RANK;i++)
if (x->v[i]) if (x->v[i])
return 0; return 0;
}
return 1; return 1;
} }
@ -189,39 +188,26 @@ static inline void bv_set_range(bitvec_t *x, int start, int end)
// find least significant 0 bit starting from position start (included) // find least significant 0 bit starting from position start (included)
static inline int bv_next_zero(const bitvec_t *x, int start) static inline int bv_next_zero(const bitvec_t *x, int start)
{ {
int position, i; int position;
position = ffsll(~(x->v[start/64] | FIRSTBITS(start % 64))); position = ffsll(~(x->v[start/64] | FIRSTBITS(start % 64)));
if(position) // found zero in same chunk? if(position)
return (start/64)%64 + position - 1; return (start/64)*64 + position - 1; // found zero in same chunk
for(i = start/64 + 1, position = 0; i < BV_QWORD_RANK && position == 0; i++) for(int i = start/64 + 1; i < BV_QWORD_RANK; i++) {
position = ffsll(~x->v[i]); position = ffsll(~x->v[i]);
if(position) // found a 0
return i*64 + position - 1;
}
return position ? (i - 1)*64 + position - 1 : BV_QWORD_RANK; return BV_QWORD_RANK; // found nothing
} }
/* static inline int bv_lowest_bit_set(bitvec_t x) */ static inline void bv_copy(const bitvec_t *from, bitvec_t *to)
/* { */ {
/* int i=0, k; */ for(int i = 0; i < BV_QWORD_RANK; i++)
/* for (i=0;i<BV_QWORD_RANK;i++) { */ to->v[i] = from->v[i];
/* if (k=ffsll((long long)x.v[i])) */ }
/* return 64*i + k - 1; */
/* } */
/* return -1; */
/* } */
/* static inline int bv_lowest_bit_set(bitvec_t x) */
/* { */
/* int k; */
/* for (k=0;k<BV_BIT_RANK;k++) { */
/* if (bv_get_bit(x,k)) */
/* return k; */
/* } */
/* return -1; */
/* } */
#endif /* __BITVEC_H__ */ #endif /* __BITVEC_H__ */

224
thickenings.c
View File

@ -578,93 +578,6 @@ typedef struct {
bitvec_t *principal_neg; bitvec_t *principal_neg;
} enumeration_info_t; } enumeration_info_t;
// calculate transitive closure; that is, fill current_level in every spot which must be marked with the current head (but was not already marked before), and -current_level in every opposite spot (including opposite to the head)
static int transitive_closure(const enumeration_info_t *info, signed char *level, int head, int current_level)
{
int is_slim = 1;
queue_t queue;
int current;
edgelist_t *edge;
queue_init(&queue);
level[info->graph[head].opposite] = -current_level;
queue_put(&queue, head);
for(int i = head + 1; level[i] != HEAD_MARKER && i < info->size; 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[i] == current_level) {
is_slim = 0;
break;
} if(level[i] == 0) {
level[i] = -current_level;
level[info->graph[i].opposite] = current_level;
queue_put(&queue, info->graph[i].opposite);
}
}
if(is_slim) {
while((current = queue_get(&queue)) != -1) {
edge = info->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[info->graph[edge->to].opposite] = -current_level;
queue_put(&queue, edge->to);
}
edge = edge->next;
}
}
}
return is_slim;
}
static inline void output_thickening(const enumeration_info_t *info, signed char *level, int current_level, int is_slim, int is_fat, long count)
{
// if printstep is set accordingly, write state to stderr
if(is_slim && is_fat && info->printstep > 0 && (count + 1) % info->printstep == 0) {
print_thickening(info->rank, info->size, level, current_level, info->alphabet, stderr);
fprintf(stderr, "\n");
}
else if(info->printstep < 0) {
print_thickening(info->rank, info->size, level, current_level - !is_slim, info->alphabet, stderr);
fprintf(stderr, " ");
if(is_slim) {
fprintf(stderr, "S");
if(is_fat)
fprintf(stderr, "F");
}
fprintf(stderr, "\n");
}
}
/* new algorithm description:
arguments: pos, neg, head (all readonly)
- newpos = union (principal_pos[head], pos)
- newneg = union (principal_neg[head], neg)
- intersection(newpos, newneg) == 0 ?
- no:
- not slim, return
- yes:
- unknown = neg(union(newpos, newneg))
- unknown == 0 ?
- yes:
- balanced, count++, return
- no:
- for all newhead in (1s of unknown):
- enumerate(newpos, newneg, newhead)
needed bitwise ops: union, intersection, negation, empty?, iteration of 1s
is there an sse op to find position of most/least significant 1?
how can we handle negatives right to head?!
*/
/* /*
ok this is screwed up, let's start over: ok this is screwed up, let's start over:
@ -693,21 +606,14 @@ static long enumerate_tree(const enumeration_info_t *info, const bitvec_t *pos,
if(next_neg != info->size/2) { if(next_neg != info->size/2) {
bv_union(&info->principal_pos[info->size - 1 - next_neg], pos, &newpos); bv_union(&info->principal_pos[info->size - 1 - next_neg], pos, &newpos);
bv_union(&info->principal_neg[info->size - 1 - next_neg], neg, &newneg); bv_union(&info->principal_neg[info->size - 1 - next_neg], neg, &newneg);
} else { // or, if there is no next_neg, just copy
bv_copy(pos, &newpos);
bv_copy(neg, &newneg);
} }
// add the range from first_unknown to next_neg to newpos // add the range from first_unknown to next_neg to newpos
bv_set_range(&newpos, first_unknown, next_neg); // including the start, excluding end bv_set_range(&newpos, first_unknown, next_neg); // including the start, excluding end
bitvec_t test1, test2;
bv_print(stderr, &newpos, info->size/2);
fprintf(stderr, " ");
bv_print(stderr, &newneg, info->size/2);
fprintf(stderr, " ");
bv_union(&newpos, &newneg, &test1);
bv_print(stderr, &test1, info->size/2);
fprintf(stderr, " %d\n", bv_disjoint(&newpos, &newneg));
// check if this leads to any conflicts (equivalently, violates slimness) // check if this leads to any conflicts (equivalently, violates slimness)
if(!bv_disjoint(&newpos, &newneg)) if(!bv_disjoint(&newpos, &newneg))
return 0; return 0;
@ -717,6 +623,7 @@ static long enumerate_tree(const enumeration_info_t *info, const bitvec_t *pos,
// do we know everything already? we have a balanced ideal then // do we know everything already? we have a balanced ideal then
if(bv_full(&known, info->size/2)) { if(bv_full(&known, info->size/2)) {
fprintf(stderr, "Found balanced ideal: "); fprintf(stderr, "Found balanced ideal: ");
bv_print(stderr, &newpos, info->size/2); bv_print(stderr, &newpos, info->size/2);
fprintf(stderr, " "); fprintf(stderr, " ");
@ -728,101 +635,20 @@ static long enumerate_tree(const enumeration_info_t *info, const bitvec_t *pos,
next_next_neg = next_neg; next_next_neg = next_neg;
while(next_next_neg < info->size/2) { while(next_next_neg < info->size/2) {
next_next_neg = bv_next_zero(&known, next_next_neg + 1); // this could return info->size/2, but that's fine for enumerate_tree int tmp = bv_next_zero(&known, next_next_neg + 1); // this could return info->size/2, but that's fine for enumerate_tree
if(tmp <= next_next_neg) {
fprintf(stderr, "%d <= %d\n", tmp, next_next_neg);
bv_print(stderr, &known, info->size/2);
fprintf(stderr, "\n");
exit(-1);
}
next_next_neg = tmp;
count += enumerate_tree(info, &newpos, &newneg, next_neg + 1, next_next_neg); count += enumerate_tree(info, &newpos, &newneg, next_neg + 1, next_next_neg);
} }
}
/*
static long enumerate_tree(const enumeration_info_t *info, const bitvec_t *pos, const bitvec_t *neg, int head, int oldhead)
{
bitvec_t newpos, newneg, known;
int newhead;
long count = 0;
if(head < info->size/2) {
fprintf(stderr, "head is %d!\n", head);
}
if(head != -1) {
bv_union(&info->principal_pos[head], pos, &newpos);
bv_union(&info->principal_neg[head], neg, &newneg);
}
for(int i = oldhead - 1; i > head && i >= info->size/2; i--)
bv_set_bit(&newpos, info->size - 1 - i);
if(!bv_disjoint(&newpos, &newneg))
return 0;
bv_union(&newpos, &newneg, &known);
if(bv_full(&known, info->size/2)) {
// found a balanced ideal
fprintf(stderr, "Found balanced ideal: ");
bv_print(stderr, &newpos, info->size/2);
fprintf(stderr, " ");
bv_print(stderr, &newneg, info->size/2);
fprintf(stderr, "\n");
return 1;
}
newhead = head;
while(newhead >= info->size/2) {
newhead = bv_last_zero_before(&known, newhead); // this will return -1 at some point, which stands for the possibility that all remaining choices are "no"
count += enumerate_tree(info, &newpos, &newneg, newhead, head);
}
return count; return count;
} }
static long enumerate_tree(const enumeration_info_t *info, signed char *level, int current_level, int head)
{
ERROR(current_level >= HEAD_MARKER, "HEAD_MARKER too small!\n");
level[head] = HEAD_MARKER;
int is_slim = transitive_closure(info, level, head, current_level);
int is_balanced = 0;
int count = 0;
// we have a candidate, check if it is a balanced thickening; if so, write to stdout
if(is_slim) {
is_balanced = 1;
for(int i = head - 1; i >= 0; i--)
if(level[i] == 0)
is_balanced = 0;
}
// comment this out (or just put it inside the if block) to save 1/3 of the runtime
output_thickening(info, level, current_level, is_slim, is_balanced, count);
if(is_slim) {
if(is_balanced) {
count++;
fwrite(level, sizeof(signed char), info->size, info->outfile);
} else {
for(int i = head - 1; i >= 0; i--)
if(level[i] == 0)
count += enumerate_tree(info, level, current_level + 1, i);
}
}
// clean up
level[head] = 0;
for(int i = 0; i < info->size; i++)
if(level[i] >= current_level && level[i] != HEAD_MARKER || level[i] <= -current_level)
level[i] = 0;
return count;
}
*/
long enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, int size, const char *alphabet, FILE *outfile) long enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, int size, const char *alphabet, FILE *outfile)
{ {
signed char *level; signed char *level;
@ -849,7 +675,8 @@ long enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, int s
if(graph[i].opposite == i) if(graph[i].opposite == i)
return 0; return 0;
if(info.size > 64*BV_QWORD_RANK) // we can only handle bitvectors up to 64*BV_QWORD_RANK bits
if(info.size > 128*BV_QWORD_RANK)
return -1; return -1;
// generate principal ideals, needed bitvec operations: bv_clear, bv_set_bit, bv_get_bit // generate principal ideals, needed bitvec operations: bv_clear, bv_set_bit, bv_get_bit
@ -882,29 +709,6 @@ long enumerate_balanced_thickenings(semisimple_type_t type, node_t *graph, int s
info.principal_pos = principal_pos; info.principal_pos = principal_pos;
info.principal_neg = principal_neg; info.principal_neg = principal_neg;
bitvec_t tmp;
for(int i = 0; i < info.size; i++) {
fprintf(stderr, "Principal ideal %2d: ", i);
bv_print(stderr, &principal_pos[i], info.size/2);
fprintf(stderr, " ");
bv_print(stderr, &principal_neg[i], info.size/2);
fprintf(stderr, " ");
bv_intersection(&principal_pos[i], &principal_neg[i], &tmp);
bv_print(stderr, &tmp, info.size/2);
fprintf(stderr, "\n");
}
// test these functions a little
bitvec_t test;
bv_clear(&test);
bv_set_range(&test, 0, 4);
bv_print(stderr, &test, 200);
fprintf(stderr, " %d\n", bv_full(&test, 6));
// don't do enumeration
// return count;
// enumerate balanced ideals // enumerate balanced ideals
bitvec_t pos, neg; bitvec_t pos, neg;
bv_clear(&pos); bv_clear(&pos);